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Abstract 

The multipolar-post-Minkowskian approach to gravitational radiation is applied to the problem 
of the generation of waves by the compact binary inspiral. We investigate specifically the third 
post-Newtonian (3PN) approximation in the total energy flux. The new results are the computation 
of the mass quadrupole moment of the binary to the 3PN order, and the current quadrupole and 
mass octupole to the 2PN order. Wave tails and tails of tails in the far zone are included up to 
the 3.5PN order. The recently derived 3PN equations of binary motion are used to compute the 
time-derivatives of the moments. We find perfect agreement to the 3.5PN order with perturbation 
calculations of black holes in the test-mass limit for one body. Technical inputs in our computation 
include a model of point particles for describing the compact objects, and the Hadamard self-field 
regularization. Because of a physical incompleteness of the Hadamard regularization at the 3PN 
order, the energy flux depends on one unknown physical parameter, which is a combination of 
a parameter A in the equations of motion, and a new parameter 9 coming from the quadrupole 
moment. 
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I. INTRODUCTION 



Inspiralling compact binaries are systems of two neutron stars and/or black holes under- 
going an adiabatic orbital decay by gravitational radiation emission. These systems consti- 
tute an important target for the gravitational-wave detectors like LIGO and VIRGO. The 
currently favoured theory for describing the binary inspiral is the post-Newtonian approxi- 
mation. Since inspiralling compact binaries are very relativistic the Newtonian description 
(corresponding to the quadrupole approximation) is grossly inadequate for constructing the 
theoretical templates to be used in the signal analysis of detectors. In fact, from several 
measurement-accuracy analyses it follows that the third post- 

Newtonian (3PN) approximation, corresponding to the order 1/c 6 when the speed of light 
c — > +00, constitutes a necessary achievement in this field. Note that the 3PN approxima- 
tion is needed to compute the time evolution of the binary's orbital phase, that depends via 
an energy balance equation on the total gravitational- wave energy flux. The energy flux is 
therefore a crucial quantity to predict. 

Following earliest computations at the 1PN level 12 1 (at a time where post-Newtonian 
corrections were of purely academic interest Vthe e nerg y flux generated by compact binaries 
was determined to the 2PN order 1 131 Il4l . Il5l llal . by means of a formalism based on 

20. Ell , and independently 



313,13 



(see also Refs. 



24, 



multipolar and post-Minkowskian approximations |l7l . Il8l [19 
using a direct integration of the relaxed Einstein equations 

EB|)- Since then the calculations have been extended to include the non-linear effects of tails 
at higher pos t-Newtonian orders. The tails at the 2.5PN and 3.5PN orders were computed 
in Refs . hd . I27I ] (this extended the computation of tails at the dominant 1.5PN order 
I2H I29L l3fll|). and the contribution of tails generated by the tails themselves (so-called "tails 
of tails") at the 3PN order were obtained in Ref. j^. However, unlike the 1.5PN, 2.5PN and 
3.5PN orders that are entirely composed of tail terms, the 3PN approximation involves also, 
besides the tails of tails, many non-tail contributions coming from the relativistic corrections 
in the multipole moments of the binary. 

The present paper is devoted to the computation of the multipole moments, chiefly the 
quadrupole moment at the 3PN order, in the case where the binary's orbit is circular (the 
relevant case for most inspiralling binaries). We reduce some general expressions for the 



multipole moments of a slowly-moving extended system [21| to the case of a point-particle 
binary at the 3PN order. The self-field of point-particles is systematically regularized by 
means of Hadamard's concept of "partie finie" 0,0, The time-derivatives of the 3PN 
quadrupole moment are computed with the help of the equations of binary motion at the 3PN 
order in harmonic coordinates (the coordinate system chosen for this computation). The 
3PN equations of motion have been derived recently by two groups working independe ntly 
with different methods: ADM-Hamiltonian formulation of general relativity 

and direct post-Newtonian iteration of the field equations in harmonic coordinates 
There is complete physical equivalence between the results given by the 
two approaches [38|,|44j. We shall find that our end result for the energy flux at the 3.5PN 
order is in perfect agreement, in the test-body limit for one body, with the result of black- 
hole perturbation theory, which is currently known up to the higher 5.5PN approximation 
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45, 4(| 4?| (see Ref. B] for a review). In a separate work Bl we report the computation of 
the 3.5PN-accurate orbital phase which constitutes the crucial component of the theoretical 
template of inspiralling binaries. 

One conclusion of the investigation of the equations of motion of compact binaries is that 
from the 3PN order the model of point-particles (described by Dirac distributions) might 
become physically incomplete, in the sense that the equations involve one undetermined 
coefficient, astatic in the ADM-Hamiltonian formalis m B- 13II I3II 37, 38| (see, however, B l) 
and A in the harmonic-coordinate approach B> 41 , kd I43I H) • Technically this is due to 
some subtle features of the self-field regularization a la Hadamard. In the present paper, 
we shall be lead to introduce a second undetermined coefficient, called 9, coming from our 
computation of the 3PN quadrupole moment. However, we shall find that the total energy 
flux contains only one unknown parameter, which is a certain linear combination of 9 and A 
entering the 3PN coefficient. All other terms in the flux up to the 3.5PN order are completely 
specified. 

The plan of this paper is as follows. Sections II to IV are devoted to the basic expressions 
of the moments we shall apply. Section V presents the needed information concerning our 
point-particle model, and Sections VI to IX deal with the computation of all the different 
types of terms in the required multipole moments. Section X explains our introduction of 
the ^-ambiguity. Finally we present our results for the moments and energy flux in Sections 
XI and XII. The intermediate values for all the terms composing the moments in the case 
of circular orbits are relegated to Appendix A. 



II. EXPRESSIONS OF THE MULTIPOLE MOMENTS 

In this section we give a short summary on the expressions of multipole moments in the 
post-Newtonian approximation. The moments describe some general isolated sources that 
are weakly self-gravitating and slowly-moving, i.e. whose internal velocities are much smaller 
than the speed of light: ?; <c. In this paper we order all expressions according to the formal 
order in 1/c, and we pose 0(n) = 0(l/c n ). In addition, the moments are a priorivalid only 
in the case where the source is continuous (for instance a hydrodynamical fluid); however, 
we shall apply these moments to the case of point-particles by supplementing the above 
expressions with a certain regularization ansatz based on Hadamard's concept of "partie 
finie" BBS. We adopt a system of harmonic coordinates, which means 



d v bT = , (2.1a) 
hT = \g\ 1/2 g^ -V^ , (2.1b) 

where g^ u and g denote respectively the inverse and the determinant of the covariant metric 
g^v, and where rf v denotes the Minkowski metric with signature +2. The Einstein field 
equations, relaxed by the harmonic-coordinate condition, take the form of d'Alembertian 
equations for all the components of the field variable, 
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UW V = -^V" , (2.2a) 



c 4 



c 4 



= \g\T^ + -A"" , (2.2b) 

167rG 

where □ = r]^ u d p d v and where we have introduced the effective stress-energy (pseudo) tensor 
t^ v of the matter and gravitational fields in harmonic coordinates. The matter stress-energy 
is described by and the gravitational stress-energy by the non-linear interaction term 
A Miy . The latter is given in terms of the metric by the exact expression 

A"" = - brti^br + d p h^d a h^ + \g^g pa d x h^d T \f x 

- g^gardxhTdph^ - g vp g aT d x h» T d p h° x + g pa g XT d x h^d T h^ 

+ \{2gWg™ - g^gn^gxrg,, - gTe9x*)d p h x *d ff h n . (2.3) 

Both the matter and gravitational contributions in r pv depend on the field h, with the 
gravitational term A Ml/ being at least quadratic in h and its space-time derivatives. 

The multipole moments of slowly-moving sources are in the form of some functionals 
of the (formal) post-Newtonian expansion of the pseudo-tensor r M ^; we denote the formal 
post-Newtonian expansion with an overbar, so t^ u = PN(r^). It is convenient to introduce 
the auxiliary notation 



= — ; Eij = . (2.4) 



C l c 



From a general study [^J of the matching between the exterior gravitational field of 
the source and the inner post-Newtonian field, we obtain some "natural" definitions for 
the Zth order mass-type (II) and current-type (Jl) multipole moments of the source. The 
physics of the isolated source, as seen in its exterior, is extracted from these multipole 
moments when they are connected, in a consistent way, to the observables of the radiative 
field at (Minkowskian) future null infinity, given in this formalism by the so-called radiative 
multipole moments. The connection between I L and Jl and the mass-type (Ul) and current- 
type (Vl) radiative moments at infinity involves up to say the 3.5PN order many tail effects 
and even a particular "tail-of-tail" effect arising specifically at 3PN. All these effects are 
known j2?| and therefore will not be investigated here but simply added at the end of 
our computation in Section XII. Here we focus our attention on the reduction to point- 
particle binaries of the general source multipole moments (in symmetric-tracefree form), 
whose complete expressions are given by 

IUt) = FP/A W 'jT,fe{ft W faS- ^[^{rf U^at, 

2(2/ + 1) 



c 4 (Z + l)(Z + 2)(2/ + 5) 



n + 2\z)xi jL ^ij } (x,t + z|x|/c) , (2.5a) 
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Mt) = FP e ab<il J d 3 x \±\ B j dzj < y,(z)z £ _ 1>a E l 



2/ + 1 

-^+i(2;)x L _i >ac E fec [>(x,t + z|x|/c) . (2.5b) 



c 2 (/ + 2)(2/ + 3) 

Our notation is as follows. L = • • •%{ is a multi-index composed of I indices; a product 
of I spatial vectors x % = Xi is denoted xl = x ix Xi 2 ■ ■ ■ the symmetric-tracefree (STF) 
part of that product is denoted using a hat: xl = STF(x L ), for instance x^ = XiXj — ^Sy, 
Xijk = XiXjXk — \{xi$jk + XjSki + XkSij); the STF projection is also denoted using brackets 
surrounding the indices, e.g. Xij = £<jfj> = \{xiVj+XjVi) — \5ijXkVk] Sijk denotes the 

usual Levi-Civita symbol (£ooo = +1); the dots refer to the time differentiation. The matter 
densities E, £< and Ejj in (|2.5jl are evaluated at the position x and at time t + z\x\/c. The 
function Si(z) is given by 

si{z) = ^[jpa - z2 y ; £ = i • (2.6) 

This function tends to the Dirac distribution when I — > +oo. Each of the terms composing 
II and Jl is to be understood in the sense of post-Newtonian expansion, and computed 
using the (infinite) post-Newtonian series 



f 1 dz8i(z)S(x,t + z\x\/c) = .,/ 2 / + 1)!! ,M 2j (4A 3 S(x,t) . (2.7) 
J_ x v ; v 1 17 ; ^2^'!(2/ + 2j + l)!!' 1 \cdt) v ' ; v ; 

Finally the symbol FP j b = o in front of the integrals in ()2.5jl refers to a specific finite part 
operation defined by analytic continuation (see Ref. [2l| for the details). Such a finite part 
is crucial because the integrals have a non-compact support due to the gravitational contri- 
bution in the pseudo-tensor, and would be otherwise divergent at infinity (when |x| — > +oo). 
The integral involves the regularization factor |x| B = |x/r | B , where B is a complex number 
and r denotes an arbitrary length scale. It is defined by complex analytic continuation for 
any B 6 C except at isolated poles in Z, including in general the value of interest B = 0. 
We expand the integral as a Laurent expansion when B — ► and pick up the finite part (in 
short FP^ =0 ), or coefficient of the zeroth power of B in that expansion. This finite part is 
in fact equivalent to the Hadamard partie finie [3l| . 

Thus, the moments depend a priori on the constant r introduced in this analytic contin- 
uation process. This constant can be thought of as due to the "regularization" of the field 
at infinity; the moments will depend explicitly on ro when the integral develops a polar part 
at B — due to the behaviour of the integrand when |x| — > +oo. As we shall see the source 
moments start depending explicitly on ro from the 3PN order. However, we know that the 
metric is actually independent of ro (more precisely, ro cancels out between the two terms 
of the multipole expansion given by Eq. (3.11) in Ref. (2l[). Indeed, as a good check of the 
calculation, we shall see that because of non-linear tail effects in the wave zone the constant 
r is cancelled out, so the physical energy flux does not depend on it. 

To the 1PN order the expressions ()2.5al) and (j2.5b|) are equivalent to some alternative 
forms obtained earlier in Refs. and [l8j|. respectively. The multipole moments in the 
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form (|2.5|) were derived in [20] up to the 2PN order, and shown subsequently in [2l| to be 
in fact valid up to any post-Newtonian order (formal ly) . On the other hand, both ()2.5a|) 
and (j2.5b|) reduce to the expressions obtained in Ref. |50j in the limit of linearized gravity, 
where we can replace t^ u by the compact-support matter tensor T^ v (hence there is no need 
in this limit to consider a finite part). Note that the source multipole moments II and Jl 
parametrize, by definition, the linearized approximation to the vacuum metric outside the 
source j2lj], but take into account all the non-linearities due to the inner (near- zone) field 
of the source. The non-linearities in the exterior field can be obtained by some specific 
post-Minkowskian algorithm (see Ref. [2jJ for proof and details). The inclusion of these 
non-linearities permits one to relate the source moments II and Jl to the radiative ones U l 
and Vl- Some other source moments Wl, Xl, Yl and Zl should also be taken into account 
(see Ref. 



21 



for discussion), but these parametrize a (linearized) gauge transformation and 
do not contribute to the radiation field up to a high post-Newtonian order. We shall check 
that these moments do not affect the present calculation. 



III. DEFINITIONS OF POTENTIALS 



Our first task is to work out the expressions ()2.5j) to the 3PN order in the case of II 
and 2PN order in the case of Jl- In this paper we shall use some convenient retarded 
potentials, and then, from these, the corresponding "instantaneous" potentials. For insertion 
into the pseudo-tensor t^ v (and, most importantly, its gravitational part A^) we need the 
components of the metric h developed to post-Newtonian order (9(8,7,8). By this we 
mean 0(8) = 0(l/c 8 ) in the 00 and ij components of the metric, and 0(7) = 0(l/c 7 ) in 
the Oi components. With this precision the metric reads 



h°° = -\v - \ (W + AV 2 ) - ^ ( Z + 2X + VW + ^V 3 ] + 0(8) , (3.1a) 
h° l = - -I (ki + Wl) + 0(7) , (3.1b) 



c 3 c 5 



F = -1 [Wit - ^W) - H (4 - l 5ij Zj + 0(8) . (3.1c) 

The potentials are generated by the components of the matter tensor T^ u or, rather, using 
a notation similar to ()2.4|) . by 

a = ; a t = — ; = T ij . (3.2) 

c 2 c 

The potential V is a retarded version of the Newtonian potential and is defined by the 
retarded integral D^ 1 acting on the source a, 



V/(x, t) = D^i-AnGa} = G J ^^^(y, t - |x - y\/c) . (3.3) 

To the 1PN order we have the potentials Vi and (together with the spatial trace W = 
Wu), which are generated by the current and stress <jj and respectively, 
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Vi = □ fi 1 {-4 7 rG'a l } , (3.4a) 

% = □ fi 1 {-47rG(^ - 8 i3 a kk ) - diVdjV} . (3.4b) 

To the 2PN order there are the potentials Ri, Z i3 , X (and also Z = Zu) whose expressions 
read 

X = n£{ - AnGa u V + W l3 d 2 3 V + Wfim 

+ Vd 2 t V + ^(d t V) 2 - 2diVjdjV^ , (3.5a) 

Ri = - AnGiaiV - aVi) - 2d k Vd % V k - ^d t Vd l V^ , (3.5b) 

Z^ = - AnG(a i3 - 5 i3 a kk )V - 2d {i Vd t V j} 

+ d t V k d 3 V k + dkVAVj - 2d (i V k d k V 3) 

- S i3 d k V m (d k V m - d m V k ) - ^{dtV) 2 ] . (3.5c) 

Next we expand the retardations and define some associated instantaneous potentials. The 
highest-order expansion is needed for the ^-potential, up to 0(5), while 0(3) is sufficient 
for Vi and Wi 3 . We write these expansions in the form 

^ = ^^-§^ + ^ + 0(5), (3.6a) 

V t = Ui + ^d 2 Xi + 0(3) , (3.6b) 

% = U l3 - |^ + ±d 2 tXl3 - \ 2 K l3 + (9(3) , (3.6c) 
where the instantaneous potentials are given by the Poisson-type integrals 

U = A- 1 {-47rGa} = G J ^l— a {y,t) , (3.7a) 

Ui = A-'i-AnGa,} , (3.7b) 

Uij = A- 1 {-47rG((T iJ - 5 i3 <7 kk ) - d t Ud 3 U} , (3.7c) 

x = 2A- 1 U = G J d 3 y |x - y\a(y, t) , (3.7d) 

X, = 2A~ 1 U t , (3.7e) 

Xij = 2A~ 1 U ij , (3.7f) 

P = 24A~ 2 U = G [ d 3 y |x - y|'V(y, t) , (3.7g) 



^ = A- l {d {i Ud 3) d 2 X } ■ (3.7h) 
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In addition, the Newtonian precision 0(1) is required for the other potentials X, B4 and Z^. 
For simplicity in the notation, we shall keep the same names for the Newtonian approxima- 
tions to these potentials, henceforth re-defined as 



Ri = A 



X = A~ 



A 



1 



1 



{ - AnGcruU + U l3 d%U + 2UidAU 

+ Ud 2 t U + ^(d t Uf - 2d i U j d j U i } , 
{ - AnGi^U - aUt) - 2d k Ud l U k - ^UdiU} 

I - 4irG(<Tij - 5ij(T kk )U -2d(iUd t Ufi 
+ diU k djU k + d k Uid k Uj - 2d(iU k d k Uj) 

- Sijd k U m (d k U m - d m U k ) - -5ij(d t U) 2 } . 



(3.8b) 



(3.8a) 



(3.8c) 



Finally the "odd" terms in (|3.6J) (having an odd power of 1/c in factor) are simple functions 
of time parametrized by 



(Beware that Q 7^ Qa.) 

IV. NOMENCLATURE OF TERMS 

The post-Newtonian metric (|3.1|) is inserted into the pseudo-tensor (|2.2b|) . in which no- 
tably the term A^, given by Eq. ()2.3|) . is developed up to quartic order h A . Making use 
of the formula ()2.7j) we obtain the source moments iz,(£) and Jz,(£) as some functionals of 
all the retarded potentials, and, then, of all the "instantaneous" potentials defined by ()3.6|) - 
()3.9|1 . We transform some of the terms by integration by parts, being careful to take into 
account the presence of the analytic continuation factor \x.\ B . The surface terms are always 
zero by analytic continuation (starting from the case where the real part of B is a large 
negative number). Notice that we use the Leibniz rule, which is surely valid in the case of 
potentials corresponding to smooth ( "fluid" ) sources. However, when we shall insert for the 
potentials some singular expressions corresponding to point-like particles, and shall replace 
the derivatives by some appropriate distributional derivatives, the Leibniz rule will no longer 
be satisfied in general. This will be a source of some indeterminacy discussed in Section X. 

We find that the moments are quite complicated, so it is useful to devise a good nomen- 
clature of terms. First, we distinguish in Ii and Jl the contributions which are due to the 
source densities S, Sj and S^- [see Eq. (J2.5J) ]. and we refer to them as scalar (S), vector (V) 
and tensor (T) respectively. Furthermore, we split each of these contributions according to 




(3.9b) 



(3.9a) 
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the value of the summation index j in Eq. ([2.7]) : for instance the S-type term denoted SI 
is defined by the set of terms in II coming from the "scalar" E in which we have used the 
formula ([2.7]) with only the contribution of the index j = (there are no S-type terms in 
Jl)] similarly we denote by SII, using Roman letters, the S-terms corresponding to j = 1 
(these terms involve a factor x 2 and a second time-derivative); and for instance VII denotes 
the set of terms in both II and Jl coming from the "vector" Sj and which have j = 2. With 
this notation the mass moment to the 3PN order can be written as: 



I L = SI+SII+SIII+SIV 

+VI+VII+VIII+TI+TII + 0(7) . (4.1) 

For simplicity's sake we omit writing the multi-index L on each of these separate pieces 
(there can be no confusion from the context). Second, the numerous terms are numbered 
according to their order of appearance in the following formulas. For instance the piece SI 
which is part of the mass moment ([4.1]) will be composed of the terms SI(1), SI(2), etc; 
similarly VII is made of terms VII(l) and so on. The numbering of terms is indicated in 
round brackets at the top of each term. The explicit expressions of all the separate pieces 
forming II is as follows. 



(2) (3) (4) 



si = fp / d^\x\ B x L h - ^A(t/ 2 ) + - -4numu 

b=o J 2nGc z c 4 7rGc 4 

(5) (6) (7) (8) 

(9) (10) (11) (12) 

(13) (14) (15) (16) (17) 

+ ^ a U a + io-l/ 3 + ly a dh + ±U jk a jk - J^UAdh 

(18) (19) (20) 

(21) (22) (23) (24) 

(25) (26) (27) (28) 

(29) (30) (31) (32) 

C331 (341 (351 (361 



;rC/r" 

(33) (34) (35) (36) 

+ -^-diUjdtUij + -^—UUidtdiU - -L-Z^U - -^-UdiUdtUi 
irGc b ixG& ixG& J nGc b 
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(37) (38) (39) (40) 



irGc^ L 1 ' 1 nGcS 1 J J 1 37rGc 6 v ' vrGc 6 

(41) (42) (43) 

(44) (45) (46) (47) 



(48) (49) (50) . 

(2) (3) 



(i) 4£7 2 



c 



A 



(4) (5) (6) 

1 'Vt^tt 1 ,*m2 2 



7rGc 4V iJ y J 2vrGc 4V 7 ttGc 4 

(7) (8) 

- ^*lU 2 - ^m(U 2 )\*\ 2 x L - U 2 d^\ 2 x L )} 

(9) (10) 
0/ _|_ Q 1 

- ^^US^x - ^jm(U%x)\*?*L - Ud 2 xW 2 x L )\ 

(11) (12) 
97 + 3 1 

" ^G^^ " 2^ d ^ UUaa ^ 2£L ~ UU °° d *(\ x \ 2S: ^ 

(13) (14) 

l r f (1) 

SHI = —77^ — ry—: FP / d 3 X | x| B <9 f 4 <^ |x| 4 X L a 

8c 4 (2/ + 3)(2/ + 5) b=o J 11 * [' 1 

(2) (3) 

- ^^WZlV 2 - ^m(U 2 )\x\ 4 x L - U 2 d t {\^x L )] } , (4.2c) 

i r {1) 

SIV = - Rl , — — — - FP / d 3 x|x| B |x| 6 x L aV , (4.2d) 

48c 6 (2/ + 3)(2/ + 5)(2/ + 7) b=o J 11 11 * V ; 

c 2 (l + 1)(2/ + 3) b=o J 

( (2) (3) (4) (5) 

(i) 9 2 1 3 

x ^ + -<7 a /J - -aU a + —d k Ud a U k + ^-^d t Ud a U 

(6) ( 7 ) (8) 

1 a fT TTT ^ , ff ao2,, , 2 _ r r2 



A(C/C/ a ) + + -a a ^ 



27rGc 2 v ay c 4 t/X c 4 
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(9) (10) (11) (12) 

1 4-2 2 

l ad tXa - -^°Ra + ~jU a a ss + -jU ak a k 

c 4 c 4 c 4 c 4 

(13) (14) 

2 1 „ 

+ ~^U k a ak + 27T Gc 4 dkUdadtXk 

(15) (16) 



(17) (18) (19) (20) 

(21) (22) (23) 

(24) (25) (26) (27) 

d k Ud a R k - -±jU u d\ x u a + ^jd t U ak d k U 

(28) (29) (30) (31) 

(32) (33) (34) (35) . 

(2) (3) (4) 



"(1) 



(6) (7) . 

- 2^ + if^PT 5)(2Z + 7) S / ^ - < 4 " 2 «> 

c 4 (/ + l)(/ + 2)(2/ + 5)i=oy dX|X| 

r (2) (3) (4) 

I (i) 1 4 4 

x <9 t 2 < a ab + j-^d a Ud b U + -^o-atf/ - ^(aU b) 

(5) (6) (7) 

+ *k* d <* udb & x + db 9 ^^ - ^h d * UkdbUk 

(8) (9) v 

+ ^A« U A U » - ^^ U M\ . ( 4 - 2h ) 



n 



Til = -Sri W7 2 ^tl 77-, r FP / d 3 x\±\ B x abL \x\ 2 

c 6 (/ + l)(/ + 2)(2/ + 5)(2/ + 7) b=oJ 1 1 abL] 1 

x dfl ai + —d a Ud h U \ . (4.2i) 

In the case of the 2PN current moment we write similarly 

J L = VI+VII+VIII+TI+TII + 0(5) . (4.3) 

The expressions of these separate pieces have the same structure as the corresponding V- 
and T-terms in the 3PN mass moment II- The differences lie only in the over-all coefficient, 
in the number of time-derivatives, and in the presence of a Levi-Civita symbol. We have 

VI = -Fp/^|x|*e^*^ (4.4a) 
VII = - 2c2(2 ) + 3) J? J d 3 x\i\ B e ab<ll x L _ 1>b 9 t 2 {same as in (JO])} , (4.4b) 
VIII = - 8c 4( 2 / + 3)(2f + 5) £o / ^I^M 4 ^* £z,-i>6<9 t 4 {same as in flOi } ,(4.4c) 



TI = 2n 21 ^ tli I o^ FP / d 3 x\±\ B e ac<ll x L _ 1>bc d t \same as in (j£2n)| , (4.4d) 

Tn = m + 2)(a + + 1 3) ( 2i + 5) is / d3x I'l'w^ as in g»> e » 

We explained that we denote the terms in the previous formulas by SI(1), SI(2), 
SI(50), SII(l), TII(2). Our convention is that this notation means that the terms involve 
their complete coefficient in front; for instance, 

SI(5) = FP / d 3 x\5c\ B x L UiMU , (4.5a) 

7rGc 4 b=o J 

sn(i4) = - 3 „ G J 2l + 3) jpJ<e*w B %a i [a i (u*)\x\ 1 x L -u 3 a i (\x\'x L )] , (4.5b) 

TI(7 > = - ,G<*(i TWTWi + 5) Ej WWgfrWW/.} ■ (4.5c) 

The notation means also that the terms include all the post-Newtonian corrections relevant 
to obtain the 3PN order in the energy flux. Consistently with that order we shall have to 
compute the mass quadrupole moment 7^ to the 3PN order, the mass octupole Ujk and 
current quadrupole Jij to the 2PN order only. Look for instance at the term SI (5) given by 
Eq. (|4.5a|) : this is a 2PN term since it carries a factor 1/c 4 . Thus, in the mass quadrupole 
Iij we need to compute SI(5) with 1PN relative precision, while in the mass octupole 1^ 
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the Newtonian precision is sufficient [the term SI(5) does not exist in the current moments]. 
Note also that a term such as SII(14) given by (j4.5b|) includes in fact two terms (which come 
from an operation by parts). Furthermore, since the different pieces (of types V and T) 
composing the current moments have exactly the same structure as in the mass moments, 
we employ the same notation for these terms in both II and Jl- For instance TI(7) denotes 
both the term in the mass moment as given by Eq. (J4.5c|) and the corresponding term in 
the current moment (with a little experience there can be no confusion). Finally, in some 
cases we split the term into subterms according to the nature of a potential therein, either 
"compact" or "non-compact" potential. The compact (respectively non-compact) part of a 
potential is that part which is generated by a source with compact (non-compact) support. 
For instance the term SI(5), which contains the potential U^ given by Eq. (|3.7cjl . is naturally 
split into the two contributions 

SI(5) = SI(5C) + SI(5NC) , (4.6) 
where is replaced by its compact (C) or non-compact (NC) parts given by 



U^ + U^\ (4.7a) 
A- 1 {-4vrG((T ij - Sijakk)} , (4.7b) 
A^i-diUdjU} . (4.7c) 

We shall split similarly all the terms containing the potentials U^, x%j> Rii %ij an d X. This 
splitting into C and NC parts is fairly obvious from the expressions of the potentials: for 
instance, 

£(nc) = A -i{_ 2dkUdiUk _ IdtUdiU} . (4.8) 

When computing the terms in the moments (|4.1j) - ()4.4jl we shall separate them into various 
categories, according to the way their computation is performed. This entails introducing 
some new terminology for the various classes. For instance we shall consider the compact- 
support terms like SI(1), or so-called Y-terms made of the quadratic product of two [/-type 
potentials [examples are VI(4) and also SI(5C)], or so-called non-compact terms like SI(5NC) 
or SII(4NC). These categories of terms are defined when we tackle their computation. The 
resulting nomenclature is complicated but turned out to be useful during the explicit com- 
putation and the many associated checks, since it delineates clearly the different problems 
posed by the different categories of terms. 



U 



U ic) 

(NC) 



V. APPLICATION TO POINT-PARTICLES 

Our aim is to compute the multipole moments for a system of two point-like particles. 
One is not allowed a priori to use the expressions ()2.5|) as they have been obtained in Ref. 
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under the assumption of a continuous (smooth) source. Applying them to a system 
of point-particles, we find that the integrals are divergent at the location of the particles, 
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i.e. when x — > yi(t) or y2{t), where yi(i) and y 2 (£) denote the two trajectories. Therefore 
we must supplement the computation by a prescription for how to remove the infinite part 
of these integrals. In this paper, we systematically employ the Hadamard regularization 
31, (see Ref. J33| for an entry to the mathematical literature). The usefulness of this 
regularization for problems involving point-particles in general relativity has been shown 
by numerous works (see e.g. jH3]). Recently the properties of the Hadamard regularization 
have been re-visited and a new set of generalized functions (distributional forms) associated 
with this regularization were introduced 0, • 

The functions F(x) we need to deal with are smooth on M 3 excised of the two points yi 
and y2, and admit when r\ = |x — yi| — > (and similarly when r 2 = |x — y2 1 — » 0) a singular 
expansion of the type 



V»6N, F(x)= r ?fW+°W). (5- 1 ) 

ao<a<n 

where the coefficients if a of the various powers of r\ in the expansion depend on the unit 
direction ni = (x — yi)/r\. The powers a of T\ are real, range in discrete steps (i.e. a belongs 
to some countable set (ai)igN) and are bounded from below (ao < a)- The functions like F 
are said to belong to the class of functions T (see Ref. for precise definitions). If F and 
G belong to T so do the ordinary (pointwise) product FG and the ordinary gradient diF. 
The Hadamard "partie finie" of F at the location of particle 1 is defined as 

(F)i = J ^ foM , (5.2) 

where dfli = dfl(ni) is the solid angle element centered on y x and of direction ni. On the 
other hand, the Hadamard partie finie (Pf ) of the integral J d 3 x F, divergent because of the 
two singular points yi and y 2 , is defined by 



Pf ui U2 d 3 xF = lim { d 3 xF 

J U^O I J r 1 >u 

+ 471 y — ( -) +4nln(-) (M+l«2| .(5.3) 

The first term represents the integral on M 3 excluding two spherical volumes of radius u 
surrounding the singularities. The other terms are such that they cancel out the divergent 
part of the latter integral when u — > (the symbol 1 <-> 2 means the terms obtained by 
exchanging the labels 1 and 2). Notice the presence of a logarithmic term, which depends 
on an arbitrary constant Ui, and similarly w 2 for the other singularity. In this paper we 
shall keep the constants u\ and 112 all the way through our calculation. We assume nothing 
about these constants, for instance they are different a priori from similar constants s\ and 
S2 introduced in the equations of motion (Section II in |43| ) . We shall see that the multipole 
moments do depend on U\ and u<i (as well as on r ) at the 3PN order. 

The strategy we adopt in this paper is to insert into the source multipole moments ()2.5|) 
the following expression of the matter stress-energy tensor T^ v for two point-masses, 
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^o-nt- parti cie = ™i« (§:) ( ) *(* - yi) + 1 - 2 , (5.4a) 



i 



1 

where mi is the (Schwarzschild) mass, yi(t) the trajectory, and Vi(t) = dyi/dt the velocity 
of body 1 [with = (c, Vi)]. This stress-energy tensor constitutes a "naive" model to 
describe the particles, since the factors of the Dirac distribution have been evaluated at 
the point 1 by means of the regularization defined by Eq. (|5.2|) . However, because of the 
so-called non-distributivity of the Hadamard partie finie, other tensors are possible as well. 
In particular, we discuss in Section X the effect of choosing another stress-energy tensor, 
which is particularly natural within the context of the Hadamard regularization, and that we 
proposed in Ref. j42|. After T^ int _ particle is substituted inside them, the moments comprise 
of many divergent integrals and we define each of these integrals by means of the Hadamard 
partie finie (|5.3jl . Therefore our ansatz for applying the general "fluid" formalism to the 
ill-defined case of point-particles is 

(II) point -particle = Pf j II [^point-particle] } ' (5.5a) 
(Jl) point -particle = Pf j Jl P^point- particle] } ' (5.5b) 

where the functionals II and Jl are exactly the ones given by ()2.5|) or (j4.1|) - (j4.4j) (including 
in particular the finite part FP^ =0 at infinity). In what follows we shall carefully apply this 
prescription, but in order to reduce clutter we generally omit writing the partie-finie symbol 
Pf. 

The relative position and velocity of two bodies in harmonic coordinates are denoted by 

dx % 

= V\ - Vl ; and v l = — = v\ - v\ . (5.6) 

To the 2PN order (only needed in this paper) the relation between the absolute trajectories 
in a center-of-mass frame and the relative ones reads, in the case of a circular orbit (see e.g. 
Ref. Q), as 

y\ = x l + 0(5) , (5.7a) 

m 

—mi +3i"y 2 5m ■ . , 

y\ = — x l + 0(5) . (5.7b) 

m 

Here mi and m 2 are the two masses, with m = mi + m^, v = mxm^l m 2 (such that 
< v < 1/4) and 5m = mi — m 2 . Furthermore, 

Gm 

7 = — , 5.8 
rc 2 
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represents a small post-Newtonian parameter of order 0(2), with r = |x|, often also denoted 
ri2, the distance between the two masses in harmonic coordinates. 

When computing the multipole moments we get many terms involving accelerations and 
derivatives of accelerations. These are reduced to the consistent post-Newtonian order by 
means of the binary's equations of motion. To control the moments at the 3PN order we 
need the equations of motion at the 2PN order. For circular orbits these equations are (see 
e.g. Q) 



dw 
~dt 

u 2 



-w 2 x+(9(5) 
Gm 



1 + [-3 + v\ 7 + 



41 9 
6 + — v + v 1 
4 



7 2 + O ( 7 3 ) 



(5.9a) 
(5.9b) 



The content of these equations lies in the relation (j5.9b|) between the orbital frequency to and 
the coordinate separation r in harmonic coordinates. However, note that the precision given 
by the equations ()5.9|) is insufficient to obtain the (second and higher) time-derivatives of 
the moments at the 3PN order. Evidently for this we need the more accurate 3PN equations 
of motion. These will be given in Section XII when we compute the total energy flux [see 
Eq. ()12.3j) below]. In addition, we shall also need for some intermediate computations the 
equations of motion for general (not necessarily circular) orbits but at the 1PN order. These 
are given by 



+ 



Gm,2 



-n 



n 



2vi + A(vxv 2 ) + -{nv 2 f + 5 



+ v[4(toi) -3(ra> 2 )]} + 0(4) 



Gm% Gm 



+ 4- 



(5.10) 



(and idem for 1 <-> 2). The notation (nvi) for instance means the usual scalar product 
between the vectors n = x/r (sometimes denoted also 1112) and vi. With these preliminary 
inputs in place, we are in a position to tackle the computation of each of the terms composing 
the multipole moments (|4.1)) - (|4.4j) . 



VI. COMPACT TERMS 

In this category we consider all the terms in (J4.1|) - (J4.4j) whose integrand involves explicitly 
the matter densities a, or as a factor, and thus which extend only over the spatially 
compact support of the source. For these terms the finite part operation FPb =0 (which deals 
with the bound at infinity of the integral) can be dropped out. With the present notation 
the compact terms are 

(i) compact term at Newtonian order: SI(1); 

(ii) compact terms at 1PN order: SII(l), VI(1); 
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(iii) compacts at 2PN: SI(3), SIII(l), VI(2), VI(3), VII(l), TI(1); 

(iv) compacts at 3PN: SI(13), SI(14), SI(15), SI(16C), SII(2), SIV(l), VI(7), VI(8), VI(9), 
VI(IOC), VI(ll), VI(12C), VI(13), VII(2), VII(3), VIII(l), TI(3), TI(4), TII(l). 

As explained earlier, it is convenient, when the potential is composed of both compact and 
non-compact parts, to separate out these pieces. Thus we shall also have the compact terms 
involving the non-compact part of a potential, namely 

SI(16NC), VI(IONC), VI(12NC). 

Evidently we have to compute the "Newtonian" term SI(1) with the maximal 3PN precision, 
while for instance a term which appears at 3PN needs only the Newtonian precision. We 
devote this section to the computation of the Newtonian term SI(1), and to one example 
of a compact term with non-compact potential: SI(16NC); the computation of the other 
compact terms is similar, or does not present any difficulty, so we only list the final results 
in Appendix A. 

From the stress-energy tensor (J5.4j) we find that the matter source densities (13. 2j) are 
given by 



<r(x,t) = /2 1( J[x-yi(t)] + 1 «-> 2 , 
<Ti(x,t) = tixv[5[x - yi(t)\ + 1 2 , 
<r y (x,*) = ii X v\v{8[yi - yi(i)] + 1^2, 

where we have introduced some "effective" masses /xi and jli defined by 



(6.1a) 
(6.1b) 
(6.1c) 



p,x(t) = ni(t) 



dt\ 



i + 



(6.2a) 
(6.2b) 



These effective masses are some mere functions of time t through the dependence over the 
particle trajectories and velocities (the accelerations are order- reduced) . Notice that, had 
we used the stress-energy tensor proposed in Section V of [42| (see also the discussion in 
Section X below), we would have found that fi\ and jli depend both on time and space, as 
they contain the factor 1/ sj—g that is given at any field point x. Using the metric (|3.1j) . 
expressed in terms of the retarded potentials ()3.3j) - (j3.5j) . we find the expressions of the two 
required factors entering the effective masses ()6.2j) up to the 3PN order; namely 



-2W + 2V< 



1 



IZ-8X + AVW - 8ViVi - -V z 

o 



+ e>(8) 



(6.3a) 



dt 
dr 



1 + - 

c 2 



1 



\v 2 + \vvl - AV lV \ + ^vt 
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1 

+ ^r 



AX + AViVi - 8Riv[ + 2W ij v\v{ - l2VViv\ 



25, 



Wv[vl + -V 3 + —V 2 vl + 
6 4 



27, 



-v v l + - v > 



+ 0(8) , (6.3b) 



where the subscript 1 means that all the potentials are to be evaluated following the 
regularization (|5.2jl . In these expressions there are no problems associated with the non- 
distributivity of the Hadamard partie-fmie; that is, we can assume (FG)i = (F)i(G)i for 
this computation (see, however, Section X). Most of the regularized values of the needed po- 
tentials at 1 (for general orbits) have been computed in Ref. (see the Appendix B there). 
Here we simply report the appropriate formulas (where ryi = |yi — y2 1 , ni 2 = (yi — y 2 )/ri 2 ). 



Gm 2 
ri2 



1 

1 + ^7 



3 Gmi 
2^7 



+ 2v 2 2 - ~(n,2< + 2 



+ 



4 Gmi 
3 ri 2 c 3 



(ni2«i2) 



Gmi 



r 12 C 



11 Gmi 5 Gm 2 
2 r 12 4 r 12 



+ - 7W2) - — ^ 2 



+ ^(^12^l) 2 - ^(«12^l)(ni2U 2 ) + y(^12^2) 2 



1 

+ ^ 



2^2 - ^(^12^2) 2 ^2 + ^(^12f2) 4 

Z o 



Gm 2 
ri2 



Gm 2 



^2 + i 



-2 -+^ 2 

ri2 



Gmi 
r 12 c 2 



-n 



12 



3 1 

--(^12^1) + -(^12^ 2 ) 



+ 0(5), 

1 Gmi 4 
+ 2^ Vl 

+ 0(3), 



4r 

G 2 mim 2 



r.1.2 



12 



iJl 



+ 



' 12 

G 2 m 2 



+ -Vg 



^(ni 2 vi)ni 2 - ^(ni 2 w 2 )ni 2 



1 

-—u 



1 



' 12 
G 2 mim 2 



12 



2 t- ^(n 12 v 2 )n\ 2 

3 Gmi 1 9 

+ -vf - 

2 r 12 4 1 



+ 0(1), 

9 

2(t>if 2 ) + -f 



11 9 11 

^(ni 2 wi) 2 + -(ni 2 fi)(ni 2 w 2 ) - ^(^i2^) 2 



+ 



G 2 m 2 r 



12 



1 Gm 2 1 2 , 1/ ^ 
127^7 " 8 V2 + % {ni2V2j 



+ 0(1). 



(6.4a) 



(6.4b) 



(6.4c) 



(6.4d) 



(6.4e) 



Notice that during the computation of the potential V at the 2PN order we used the 1PN 
equations of motion for general orbits: these are given by Eq. (|5.10|) . In addition to the 
above, we need the trace W = Wu at 1PN order. [To the order considered in ()6.4cj) we have 



Wij.) By a computation similar to those of Ref. |5JJ we get 
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Gm2 
ri2 



Gm\ 1 Gm 2 
n 2 2 r i2 



-2^ 



2G 2 mirri2 



r 2 n c 



{ni2V 12 ) 



G 2 m\m2 
r 2 2 c 2 



gGmi 1 Gm 2 
ri2 2 na 
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+ o w i + -7T V 2 + ("12^1) + 2(ni 2 vi)(ni 2 v 2 ) - 2(^12^2)' 



+ 



G 2 m\ 
r 2 2 c 2 



9 ^ / \ 2 

-v 2 + -(71^2) 



+ 



Grri2 
r 12 c 2 



[-2u| + (n 12 t; 2 ) 2 t; 2 ] + 0(3) . (6.5) 



Inserting these expressions into ()6.3|) we obtain the 3PN fX\ and then straightforwardly 
compute SI(1). In the quadrupole case I — 2 it is given by 



SI(1) = J d 3 x faff'Si + 1^2 

= hyfyf + 1 <-> 2 . (6.6) 

The final result for circular orbits [using the relations (|5.7j) - (|5.8jl ] reads then 



SI(1) 



mis 

7 3 



1 + 1(1-51/) 



— (13 - 61z/ + 25z/ 2 ) 
8 



+ — (149- 573// + 354i/ 2 - 29z/ 3 ) 
16 



(6.7) 



The sensitivity of this result to the choice of stress-energy tensor for point-particles (in 
accordance with the "non-distributivity" of the partie finie) is discussed in Section X. 
Other interesting terms in this category are 



r(NC) 
' ab 



SI(16NC) = i FP y d 3 x\x\ B x L a ab U { a 
and the similar VI(IONC) and VI(12NC). Applying our computation rules we get 

4m l ah ( , , \<i, , M->r r (NC)\ 



SI(16NC) = -^vf [(y l +r l n 1 ) < \y 1 +r 1 n i y > U^ )) ) i + 



(6.8) 



(6.9) 



where we have written x l = y\ + r\n\ valid in the vicinity of the point 1. The result follows 
from applying the regularization ([5.2)1 . with the help of the Newtonian approximation of the 
NC potential. The interesting point is that the regularized factor in ()6.9|) is different from 
Ui l yi > f^ir°' ) ) as a consequence of the non-distributivity. See Section X. 



VII. QUADRATIC TERMS 

In this category we consider all the terms whose support is spatially non-compact (hence 
the finite part operation FP# =0 plays a crucial role), and which are made of the integral 
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of a product of two derivatives of compact-support potentials. Furthermore we sub-divide 
the quadratic terms into sub-categories Y-, S- and T-terms named after the functions Yl, 
Sl and Tl defined below, and we classify all these terms according to their dominant post- 
Newtonian order. The exhaustive list follows. 

(i) F-terms at 2PN: SI(4), SI(6), SI(7), SII(7), VI(4), VI(5), TI(2); 

(ii) y-terms at 3PN: SI(31), SI(35C), SI(37C), SI(38C), VI(16), VI(20), VII(6), VI(19), 
VI(21), VI(25C), VI(26C), VI(27C), VI(29C), TI(6), TI(7), TI(8); 

(iii) S-terms at 3PN: SII(3), SII(4C), SII(5), SII(6), SIII(2), VII(4), VII(5), TII(2); 

(iv) T-terms at 3PN: SI(17), SI(19C), SI(21C), SI(24), SI(25), SII(9), VI(14), VI(15), VI(17), 
VI(18), TI(5). 



The Y- and S-terms involve the product of two compact-support potentials U, Ui or U 



(c) 

ij 



while the T-terms involve a product of one of the latter potentials (of type U) and a potential 



of the type Xi Xi or Xij [ see (13. 7j) ]. Compared to F-terms, the S'-terms contain in addition 
a factor |x| 2 inside their integrand. In the two-body case these compact-support {/-type 
potentials read 

U = 9h + i^ 2 , (7.1a) 

n 

^ = ^-v\ + 1 <-> 2 , (7.1b) 

n 

Up = _ 5 ij v f) + i <_> 2 . (7.1c) 

The potentials of type x are obtained by replacing l/r\ by r\ in these expressions. Then 
from the structure ~ l/r% + l/r 2 or ~ r% + r 2 it is not difficult to express all the Y-, S- 
and T-terms with the help of three and only three types of elementary integrals Yl, Sl and 
Tl respectively (where L = ■ ■ - i\ denotes the multipolar index). Two examples in the 
quadrupole case ij are 



AG 



,a„.b 

G d 2 



SI(4) = TtnlfituM dabYij + iwl daYj] + 1 «-> 2 , (7.2a) 

C 4 2 2 



SII(4C) = ^m x m 2 ^[{vf - 6 ab v 2 2 ) d ab S tJ ] + 1 <- 2 . (7.2b) 

Since SI (4) is a 2PN term it needs the relative 1PN precision (for simplicity we do not write 
the post-Newtonian remainders). The elementary integrals are defined by 

n(yi,y 2 ) = ~ FP / d^\±\ B ^ , (7.3a) 
Itx b=o J r\Ti 

&(yi,y a ) = — FP / rf 3 x|x| B |x| 2 ^ , (7.3b) 
2% b=o J r 1 r 2 
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1 nn f j3„|~|B* r l 



T L (y 1 ,y 2 ) = -— FP / ^x|x|^^ . (7.3c) 
2n b=o J r 2 

In these definitions, the finite part at infinity is absolutely crucial (it comes directly from 



the formalism |2(j, |2l|). However, it is easily seen that the integrals are convergent near 
the two bodies so the Hadamard partie finie is not needed. The integral Yr, agrees with the 
dehnition used in 00 and is equivalent with the alternative form proposed in Ret. Q. 

We present several derivations of the closed-form expressions of these integrals for arbi- 
trary I. This permits us to introduce some techniques which are necessary when we compute 
some more complicated integrals in Sections VIII and IX. The first method consists of writing 
the multipolarity factor in the form 



p=0 



V 



where ( ) denotes the binomial coefficient (and <> refers to the STF projection). Inserting 
this into the integral Yl, it is easy to obtain the equivalent expression 

^-sSO^^^fc/^'?} • (7 ' 5) 

Next we compute the integral inside the brackets of ([7.5)1 . Let us show that the polar part of 
this integral when B — » is zero. We replace the integrand by its expansion when |x| — > oo 
(any pole at B = necessarily comes from the behaviour of the integral at infinity), we 
integrate over the angles and look for radial integrals of the type J °°<i|x| |x| B_1 which are 
the only ones to produce a pole. However these radial integrals do not exist since after the 
angular integration the powers of |x| are only of the type B + 2k where k is an integer. So 
the integral in ()7.5|) can be computed by analytic continuation down to the value B = 0. 
We obtain (Vp G N) 

FP / d 3 x = - „ 12 r » ( 7 - 6 ) 

B=oJ 11 r 2 (p + l)(2p+l) V ; 

which is a particular case of the Riesz formula , valid for any a, b E C except at some 
isolated poles: 

r (-§) r (-1) r 12 1 ; 

(T denotes the Eulerian function). A closely related reasoning to prove ()7.6|) is to replace 
the regularization factor |x| s by its expansion when B — > 0, i.e. 



|x| s = ff <j 1 + — In 



B 
~2 



1 + 2^^ + ^ 



+ 0{B 2 ) \ . (7.8) 
r 1 rfj J 

Since the integral does not develop any pole when B — > 0, the term of order B cannot 
contribute, nor any of the higher-order terms 0(B 2 ). This means that we can replace 



21 



the regularization factor \5t\ B by ff (where r± = ri/ro). From the Riesz formula, with 
a = B + 2p — 1 and b = — 1, and computation of the limit -B — > we get the same result. 
Thus, plugging ()7.6|) into (|7.5|) we find the explicit expression of Yl as 



n(yi,y2) = ri2^ ( J 



pj V + 1 



(7.9) 



where y{ 2 = y\ — y\ and r 12 = |y 12 1 - In terms of y\ and y\ the expression is simpler: 



Yr 



ri2 



<L-Q Q> 



9=0 



Using exactly the same method we find for the S^-integral, 



S, 



±Q(-)'yf l - p y?i 



p=0 

ri2 



^ ( P+ l- 2 l 

(p + 2)(p + 3)V 3 



2yi-y 



12 



(/ + !)(/ + 2) 



9=0 



<L-QyQ> 



and, for the T^-integral, 



(7.10) 



yf 



p+2 p+1 



(Z + 1 - q)yl - -{q + 1)(Z + 1 - g)y? 2 + (g + l)y 2 2 



(7.11) 



TV 



E 



1 12 

3 U \pJp+ 2 



yf L - p y^ 



12 



3(1 + !)(/ + 2) 



E(^+ 1 )^ <i " Q ^ > 



g=0 



Notice that Si, can be deduced from Tl and Yj, using the formula 



S L = {l-2y\ di)T L + ylY L . 



(7.12) 



(7.13) 



The integrals Yl, Sl and Tl vanish in the limit yi — > y 2 . As is clear from the defining 
expressions f|7.3|) there is no problem with the latter limit, in the sense that it does not 
introduce any singularity at the point 1. This justifies a posteriori our neglect of all the 
"self" contributions (proportional to m\ and m 2 ) in the quadratic terms; see the examples 
given by Eqs. ()7.2jl . However, when we compute the cubic and non-compact terms in 
Sections VIII and IX we shall find some important non-zero self contributions. 

Another method for the computation of the integrals ()7.3j) is based on the set of functions 
defined by 



g = ln(r x + r 2 + r 12 ) 



(7.14a) 
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/ = jjrO" 2 + rj- r 2 l2 ) (g-^\+ ^(r 12 ri + r 12 r 2 - r 1 r 2 ) , (7.14b) 

12 X / l\ 1 

/ = q ( r i + r i2 - r l) ( - g J + - (nr 2 + r i2 r 2 - nri 2 ) , (7.14c) 

21 1 / l\ 1 

/ = q ( r 2 + r i2 - r i) ( 9 ~ g J + g ( r i r 2 + rviTx - r 2 r 12 ) , (7.14d) 



which satisfy, in the sense of distribution theory, 



Ag = ; A l9 = ; A 2 g = , (7.15a) 

r 1 r 2 r x r X2 r 2 r 12 

Af = 2g; AJ = A 2 f = ^ , (7.15b) 

7*12 r 12 
12 r 12 12 r 

A f = -1- A X /= 2g; A 2 f=^, (7.15c) 
r 2 r 2 

21 r 21 r 21 

A / = — ; A, f= ^; A 2 /= 2g , (7.15d) 
n ri 

where the Laplacians A = <9j<9j, A x = i^i^j, A 2 = 2 d i2 di. Let us take the example of the 
integral Yl. With the help of ()7.14a|) it can be re- written as 

Y L = ~ FP / t/ 3 x |x| B x L A<7 . (7.16) 

We operate the Laplacian by parts, discard the -B-dependent surface term which is zero by 
analytic continuation, and use the formula A(|x| b £l) = B(B + 21 + 1)|x| b ~ 2 xl. Hence, 

Y L = -— FP <B(B + 21 + 1) [ d 3 x |x| B |x|~ 2 x L ^ 1 . (7.17) 
2tt b=o [J J 

Because there is an explicit factor B in front of the integral we need to look only at the polar 
part when B — > 0, which depends only on the behaviour of the integrand at the upper bound 
r = |x| — > +oo (this r should not be confused with r = T\ 2 as we sometimes denote the 
orbital separation). Thus we are allowed to replace the function g in (J7.17)) by its expansion 
at infinity. It can be checked that the (simple) pole of the integral in (|7.17j) is produced 
exclusively by the term in the expansion of g of order r^ 1 ^ 1 . Let us consider the quadrupole 
case 1 = 2. We have 



g = m(2r) + -{--. }+-{•-•} 

+ pf {"-f + (*Jfc)( n H») + (^) 2 ] +•••} + O (^j , (7.18) 

where the dots indicate some terms which yield no contribution to the present computation, 
either because they do not belong to the relevant order r~ 3 or they will be zero after angular 
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integration. Thus the formula (|7.17|) becomes in this case 



y y =Fp{-2B(fl + 5) | + °° dr r^r- 1 J ^ h/-f [(ny,) 2 + (n yi )(ny 2 ) + (ny 2 ) 2 ] } . 

(7.19) 

The notation for the radial integral means that only the bound at infinity contributes to its 
value. The latter expression is easily transformed into 

Ya = r f [yf l]> + yfy? + y? ]> ] , (7.20) 

in agreement with the more general result (|7.1()jl . The same method works for Sl as well, but 
one performs two successive integrations by parts using the functions g and /. Concerning 
Tl, one integration by parts is sufficient but using the function f 12 (labels such as 12 are 
placed at the top when the quantity appears in an equation and as right-side superscripts 
when it is within the text). 

With the latter expressions of the elementary integrals Yj,, Sl and Tl we obtain all the 
quadratic terms. The results in the case of circular orbits are displayed in Appendix A. 



VIII. CUBIC TERMS 

By cubic terms we refer to all the terms which are made of a product between three 
(derivatives of) compact-support potentials U and Ui [there are no such terms involving the 
tensor potential U^]. From (|4.2jl we can check that the only cubic terms appear at the 
3PN order. These are 

SI(26), SI(27), SI(28), SI(29), SI(30), SI(34), SI(36), SII(13),VI(22), VI(23), VI(24). 

Let us proceed in a way similar to the computation of the quadratic terms, i.e. by expressing 
the terms as functionals of some elementary integrals that are computed separately. Since 
the cubic terms are 3PN, their computation can be done using the Newtonian potentials 

U = — + 0(2) + 1^2, (8.1a) 

n 

Ui = — v\ + 0(2) + 1^2. (8.1b) 

n 

For simplicity we gather in one computation the sum of all the cubic terms in SI [and similarly 
in VI; there is only one cubic term in SII, which is SII(13)]. In the case of mass-type moments 
we get 



SI(26 + 27 + 28 + 29 + 30 + 34 + 36) 

(- 5 -°) 512 



G 2 m\ 



32 • (-3.°) 



+T v( Y L 
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G m\m<i 
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-v\v{ dij Yi 
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-Svlvidij Y L 



+ \-2v\v{ - 6v[v{ + 8S ij ( Vl v 2 ) + 8tf* J 'u?] di dj ( Yl 

_ (-4,-1) 1 

+ [15(u 1 t; 2 ) + 3v 2 } Y L ^ + 1^2, 

AG 2 d 2 f (- 3 >°) (-2,-1) 
SII(13) = — —lm\ Y L +Zm\m 2 Y L \ + 1 <-> 2 



VI(22 + 23 + 24) 



8G 2 (2/ 
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r7 f (- 5 >°) 
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3 (-2,-1) (-2,-1) 

-M-v k 2 ) dadk Y aL -vfdkdk Y aL 
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Ir rf-V%)dak Y aL 

16 i 



3 5 \ (- 4 '" 1 ) 
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(8.2a) 
(8.2b) 



(8.2c) 



In the case of the current-type moments there are only the Vl-terms, which admit a formula 
analogous to (|8.2c|) . As we see, we could express all the cubic terms by means of a single 
type of elementary integral, 



(n,p) 

Y L 



,yi,y2j 



^FP 

2% B=0 



d 3 x\±\ B x L r^r p 2 



(8.3) 



(-1,-1) 



and Tl 



of which some particular cases used in the previous section read Yl = Y 1 
Yl [ we use right-side superscripts (n,p) when the quantity appears within the text]. The 
integral ({HUH i s well-defined in the vicinity of the points yi and y 2 only when n > —3 and 
p > —3. When this is not the case - for instance the integral y[ 3 '°' ) appearing in ()8.2j) - 
one should add the Hadamard partie-finie operation Pf defined by (|5.3|) and depending a 
priori on two constants U\ and u 2 - According to our convention we generally do not write 
such parties finies, but they are always implicitly understood. 

The integral Y L ' is perfectly well-behaved near the two bodies (like Y L , Sl and T L 
considered in Section VII) , so it does not need the partie finie. We substitute in it a formula 
obtained from (|7.4jl by exchanging the labels 1 and 2, obtaining 
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FP 
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(8.4) 



Next we replace the regularization factor |x| B by its expansion around B = already written 
in (|7.8p . Since the integral can develop simple poles at most, we can limit ourselves to the 
first order in B. Then the integral in the brackets of ()8.4j) reads 
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The first term follows from the Riesz formula (J7.7|) . and the second term depends only on 
the poles developed by the integral at infinity (because of the explicit factor B in front). 
Now, contrarily to the case of the integral Yl = Y^ -1 '" 1 ) investigated in Section VII, we find 
that this second term gives a net contribution to the integral, straightforwardly obtained 
from expanding the integrand when r = |x| — > +00. The final values that we obtain in the 
quadrupole and octupole cases {1 = 2 and I = 3) of interest are 
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(8.6a) 



(8.6b) 



Note the occurence of some logarithms of ryi = ^12/^0- Applying on these values the 
point-1 Laplacian Ai = idu, and using Airf 2 = 2rf 4 (a statement valid in the sense of 
distributions), we obtain 



(-4,-1) 1 /8 ■• 4 ■ 1 

y* = 6r ( \y? 3> - \y?y? - \y? 3> ) • <' 8 --> 



12 



^ \? ^ 1/16 <i,/ c > 8 < j,- k> 2 <i ,- fc> 1 <ijfe>\ / «i \ 

Yijk = ^-5-2/1 " 5^1 ^2 " gl/f 2/2 - g2/2 J • (8.7b) 

Alternatively, the results (|8.7|) can also be obtained by the same technique as used previously 
for Y£ ' (i.e. from the Riesz formula and search for the pole at infinity). 
The computation of the integral Y^~ 3 '°\ defined by 

(-3,0) 1 r ™ r 

Yl (yi) = -=- FP / <i 3 x|x| B ^, (8.8) 

is a priori more tricky because this integral necessitates the Hadamard partie finie for curing 
the divergence at the point yi. Actually, the same method as before, based on the Riesz 
formula, could be used because we know that the Hadamard partie finie can also be obtained 
as an analytic continuation (see e.g. 41]). We prefer here to vary the techniques and to 
present some other derivations. We split the integration domain M 3 into a ball centered on 
yi with some fixed radius IZx, and the complementary domain, i.e. r\ > TZ±. The partie 
finie applies only on the "inner" domain, surrounding the singularity 1, and the finite part 
FPb=o applies only on the integral extending to infinity. Hence, 
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In the first term we recall that the partie finie depends on a constant U\ [see the definition 
(15. 3|) ]. For this term we readily find 

--Lpf Ul / rf 3 x% = -2y L x \n C*±) . (8.10) 

On the other hand, one must replace into the second term the factor \5t\ B by its 5-expansion 
as given by (|7.8|) . This yields two contributions: one is immediately computed using the 
properties of the analytic continuation, the other contains an explicit factor B and therefore 
relies on the existence of poles at infinity: 



— FP / d 3 x|x| B % 
2tt b=q J ri>ni r{ 



j3 ~B X L , 

a x r l — In 



l + 2™2 + % 
ri r\ 



(8.11) 



As expected, the sum of the two contributions ()8.10|) and (|8.1H1 is independent of the 
intermediate length scale TZ%. Indeed, the integral in the second term of ()8.11|) does not in 
fact depend on TZ\ as it depends only on the infinite bound. We obtain 



(- 3 >°) T /uA 1 



d 3 xr~f^hi 



1 + 2™! + % 



(8.12) 



The computation of the second term proceeds along the same line as for the reduction of 
Yl in (J7.17j) . We expand the log-term up for instance to the order \jr\ necessary to get the 
quadrupole case 1 = 2, 
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The integral follows immediately. This method yields the results (cases I = 2, 3) 
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The results depend on the Hadamard-regularization constant u\. 



(8.15a) 
(8.15b) 



27 



We present another derivation of the integral Y L 
distribution theory (see e.g. |33| ) 



(-3,0) 



based on the interesting formula of 
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+ 47r5(x-yi^ 



(8.16) 



[Notice the sign of the distributional term, +Att5i, opposite to the sign in the more famous 
formula A(l/ri) = —Att5i.] With Eq. f)8.16j) one can re-express 3 '°' ) in the form 
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Here the first term comes from the delta- function in (|8.16|) . Integrating the second term by 
parts, we get 
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Following the same principle as before, we compute the remaining integral by looking at 
the pole at infinity. The result is in agreement with the earlier derivation (as we checked in 

consists 

i) 



the case 1 = 2). Let us also mention that still another method to compute Y L 
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The limit is singular since Y^ 2 ' 



of taking the limit y 2 ^ ji ^ unc m^giai ± L 
diverges when the two particles merge together. In fact the limit must be taken in the sense 

we obtain the 



of the Hadamard partie finie (|5.2j) . Indeed, applying Eq. (5.5) in Ref. 
following limit relation between Y£ 3,0 ' ) and Y^ 2 ' ^: 
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Inserting for instance the result for Y^ ' obtained in ()8.6a|) we recover exactly the function 



-(-3,0) 



*3 



given by ()8.15a|) . 



is identically zero. We apply the point-1 Laplacian Ai onto the expression of Y L 
the known formula of distribution theory 



Finally it is easy to see that the function Y^ also needed into the cubic terms 
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and readily obtain, for any I, 
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. (8.21) 
The results for the cubic terms in the case of circular orbits are reported in Appendix A. 
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IX. NON-COMPACT TERMS 



The most difficult part of the present analysis is the computation of the so-called "non- 
compact" terms, which are cubically non-linear terms (like the cubic terms) made of the 
product of a compact-support potential like U and a quadratic "non-compact" potential 
like ujj* . The complete list of non-compact terms is 

SI(5NC), SI(19NC), SI(20), SI(21NC), SI(33NC), SI(35NC), SI(37NC), SI(38NC), SII(4NC), 
VI(25NC), VI(26NC), VI(27NC), VI(28NC), VI(29NC). 



A. Expressions of the NC terms 

As before, here again our strategy is to express the non-compact terms as functionals of 
certain elementary integrals, that are computed separately. We substitute inside the sources 
of non-compact terms the appropriate post-Newtonian potentials computed for two particles 
on a general orbit. The compact potentials U, Ui and (and similar expressions for the 
x's) were already given by (|7.1)1 . Here we list all the non-compact potentials needed for this 
computation [see (j3.7j) - (|3.8|) for definitions]. The potential uj^ C ^ is the only one which is 
needed at 1PN order; the other potentials are Newtonian. 
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Here, g, f, f 12 and f 21 are defined by (|7.14|) . and we denote e.g. ^ = xd^djg (see Ref. j5l| 
for the expression of igj); the acceleration is a\ = dv\/dt; the parenthesis around indices 
denotes the symmetrization (and G = 1). 

Notice that we have chosen to express the non-compact potentials by means of g, f, 
f 12 and f 21 . But these functions constitute merely some particular solutions of the Laplace 
equations (J7.15|) we have to solve, and the question arises of which solution is the correct one. 
The most general solution will be obtained by adding to the particular one an homogeneous 
term, solving a source-free Laplace-type equation. We have checked that the only possible 
homogeneous solutions, that are regular at the origin, are constants or linear functions of the 
position, and that these are always either cancelled by some spatial or time derivatives, or 
disappear at the end of our computations. This justifies our use of the particular solutions 
f!7.14|) . (Similarly, we found that the same happens in the computation of the 3PN equations 
of motion, where these particular solutions are sufficient 43].) 

The potentials ()9.1|) contain a "self" part, proportional to m 2 or m 2 (before replacement of 
the accelerations), and an "interaction" part, proportional to m 1 m 2 . Similarly the sources of 
the non-compact terms will involve a self part, proportional to mf or tttJ, and an interaction 
part, proportional to m 2 m^ or mim 2 ,. At the 2PN level, all the self parts cancelled out 
in the multipole moments [l3j. At the 3PN level, we shall find that the self parts bring 
a contribution to the moments. [Actually, we shall argue in Section X that the self parts 
are unknown.] For treating the NC terms we used the standard distributional derivative 
0,0]. Thus, we have, for instance, 



A— = -4vr5i , (9.2a) 

%(3 = ^-^. (»■*» 

A(i)^- 1 ^, (U.) 

^(ij^i^-l^-^- (9-2d) 

However, the use of the standard Schwartz derivative can be justified only when the terms 
involved are multiplied by some smooth functions. In the case of the self parts of NC terms, 
this will not be true in general, so the Schwartz derivative gives some ill-defined contributions, 
composed of the product of a delta-function and a singular function. In Section X we 
consider a well-defined way to do the computation of the self terms, which is based on 



the distributional derivatives proposed in Ref. |41[. From the discussion in Section X we 
conclude that one must add to the present computation some undetermined terms taking 
into account the ambiguities in the choice of the regularization and distributional derivatives. 
All the expressions in Eq. ()9.3)1 below are modulo these ill-defined contributions and we can 
safely proceed with the knowledge that our procedure is unambiguous and complete. We 
are securely protected from such ill-defined contributions at this stage since we shall add 
such terms with an arbitrary coefficient in Section X. We obtain the following expressions 
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of the non-compact terms, as functionals of several new types of elementary integrals (we 
pose D = \di 2<9j and G — 1). In the case of the mass-type moments: 
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We have similar expressions (involving Vl-type terms) for the current moments. The ele- 
mentary integrals parametrizing the NC terms include some generalizations of the integrals 
already introduced in Section VIII, 
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As usual the Hadamard partie finie Pf is to be added when the integral diverges near the 
particles. The logarithms in ()9.4c|) - (j9.4d|) contain the constant r through the notation 
fi = tx/tq. In addition we have the more involved integrals 



riaP 



K P l 



FP / d 3 x 

2n b=o 



FP / d 3 x 

2n b=o 



1 

2~n 

1 

2~n 

— FP / d 3 x 

2lX B=0 



— FP / ci 3 x 

27T B=0 



- — FP / d 3 x 

2n b=o 



B x L d a p | — | a g , 



x L d a p I — ) b9 , 
Tx 



B x L d P ( — ) a 



xlO p ( — ) fc^ 



X L d a P [ — ) af 



(9.5a) 
(9.5b) 
(9.5c) 
(9.5d) 
(9.5e) 



33 



12 

F[ = 

21 

F[ = 



m = 



Ql = 



-— FP / d 3 x 

2n b=o 



Iff [ 

2tt b=o J 

4- Fp / 

2n b=o J 

±fp r 

2n b=o J 



d 3 x 
d 3 x 
d 3 x 



- — FP / d 3 x 

2lt B=0 



-— FP / d 3 x 

27T _B=0 



\ \ 12 

B ^L<9 a P ( — J a/ , 

1 \ 2 l 

n 



B x L d aP ( — )„/, 



B ~ 



B x L \x\ 2 d aP ( ) a g , 



3 x L dj 



B x L d i 



d i\-\ ^9 



d i\-)i9 



(9.5f) 
(9.5g) 
(9.5h) 
(9.5i) 
(9.5j) 
(9.5k) 



The notation is e.g. a f = id a f ', g b = 2 dbg, a 9b = \d a2 d b g (notably k g k = Dg). The last two 
integrals are related to some previous ones by 
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B. Computation of the elementary integrals 

The techniques developed in Sections VII and VIII can be used to compute many of these 
integrals. Concerning S^'^ we need only the particular case I = 2 and (n,p) = (—5,0). It 
is computed by the same methods as used for Y^ -3 ' -*; we find 
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Next the group of integrals constituted by the N 

-2 



's and M^'^'s is obtained in a fashion 



similar to the one employed for 



^ in Section VIII, i.e. basically by application of the 
Riesz formula. The logarithms in these integrals are included by differentiating with respect 
to the complex parameter B. The relevant results are 
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The remaining integrals, denned by (|9.5|) . are more difficult, but we have been able to 
obtain all of them using several different methods, adapted to the computation of each of 
these integrals separately. We shall not present all the details of these computations but 
simply outline some examples. Consider the integral Kl defined by (|9.5cj) with p — 0, i.e. 



K L = ~ FP I d 3 x\5t\ B x L ^-. 

27V B=0 J T\ 

Using the fact that g/r\ is a Laplacian, 
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we can integrate by parts and transform Kl into an integral containing an explicit -B-factor, 
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From a previous argument, the value of the integral depends only on the possible occurence 
of a pole ~ 1/B at infinity. As the pole is easily computed from expanding the integrand at 
infinity, we obtain in this way the expression of Kl- Next, from the formula 
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where one should be careful about considering A\r 1 1 in the sense of distributions [i.e. 
Axrf 1 = — 47T we deduce Gl from the Laplacian of K L . Indeed, as a consequence of 
(EH, 
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and we can easily show that 2 '°' ) is actually zero. Alternatively, one can prove also that 
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This provides a check of the computation. 

To compute G S L (in the quadrupole case L = ij, say) we use a different method. We 
remark that G|- obeys a Laplace equation, with respect to the point 2, with known source: 



A 2 =da [ J dasYij . 

T\2 J 1 
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Here, Yy is known from (|7.10|1 . The right-hand-side of Eq. ()9.15|) is expanded, and we 
obtain a particular solution of this equation by integrating each of the terms. Now G\^ is 
necessarily equal to this particular solution plus some solution, regular at the origin, of the 
homogeneous equation. Taking into account the index structure of G\^ and the fact that it 
has the dimension of a length, we find that the homogeneous solution is parametrized by 
solely two numerical constants a and b. At this stage we have 



Gt, = - ^Vu j>S rli + IvlfyfrTi - ly<« <P> s lnf 12 - iy<^>' lnr~ 12 

+ ayf5 j>s + by< i 5 3>s . (9.16) 

Incidentally, this expression already gives the complete result for the gradients id s Gf, and 
2<9 s (jrf ? -, because the gradients of the homogeneous terms are zero. To compute the constants a 
and b we need some extra information, which is provided by the contracted product between 
y{ 2 and Gfj. Indeed this contraction is a known quantity thanks to the identity 

y s i2 Gtj = - U + y s 12 ds^j Gij + iAi , (9.17) 

(—2 1) 

where has just been obtained previously. Here, 1 } can be computed from the Riesz 

/ 2 _ ]\ ______ 

formula exactly like for ' in Section VIII. [When deriving ()9.17j) we take account of the 

fact that 2,0 ' ) = 0.] Comparing the result for yf 2 G|,- with the one obtained directly from 
()9.16|) we find three equations for the two unknown constants a and b. This overdetermined 
system fixes uniquely the constants to the values a = 63/100 and b = —257/900. 

The preceding method was successfully applied to several integrals of the type (j9.5j) : 
that is, we (i) compute the "source" of the Laplace equation satisfied by the integral with 
respect to the point 2 (the source is computable because A 2 applies only on the part of 
the integrand containing the functions g, f, etc., and we can make use of Eqs. (|7.15j) : 
with respect to the point 1 this would not work), (ii) compute a particular solution of this 
equation, (iii) write down the most general form of the homogeneous solution in terms of 
a few arbitrary coefficients (this works only when the dimension of the integral is a small 
power of a length so that the number of unknown coefficients is small), (iv) compute the 
coefficients using the extra information provided by the contraction with respect to yi 2 . 
Alternatively to (iv) one can use an angular average with respect to ni 2 [see (19.29)1 below]. 

As a verification let us introduce the new integral 

% = -7^ FP f rf 3 x lifxijda (-) g a . (9.18) 
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From the easily checked formula 

(A : - A) ( 2-) = 4 1- + 29 a (-) g a , (9.19) 

we deduce a relation between R^a and some computable quantities, 
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The value of the last integral comes from the pole at infinity - same method as before. 
Having obtained Rij, the verification is that 2d s G?j, which on one hand is computed from 
()9.16|) . on the other hand should be given by the following alternative expression: 
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which is obtained by some integrations by parts inside the integrand of 2d s Gfj. Of course, 
the value of computed by (|9.2U|) is such that (|9.21|) is also satisfied. 

Once Gfj is known we can deduce another needed integral, i.e. 2d s K£?, from the identity 
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+ ^-Fp|fi(fi + 5) y"rf 3 x|x| B |x|- 2 %9 a (7)2} • ( 9 - 23 ) 

Again the last integral makes no problem. Next, from both R^ and 2d s K^, we can further 
deduce 2<9a(bC° ! ). Indeed the other identity 
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Some other integrals are connected directly to the simpler y-type integrals. For instance, 
the integral (|9.5d|) is given by 
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(using the facts that Y i3 
can check that its trace 
perfect agreement with 



(-2,0) 



= Y s 



(-4,0), 



. Once the value of this integral is obtained, we 
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the final reduction being obtained thanks to the known formula (see e.g. j5l|) 

1/1 1 1 \ , 
k9k = -( • 9.28 

2 \r x r 2 nr 12 r 2 r 12 J 

Still another method is useful in our computation. All the integrals are certain functions 
of the two points yi and y 2 , and it is advantageous to consider their angular average with 
respect to the relative direction ni2 between the points, with the vector yi being fixed. As 
it turns out, the average is much easier to compute (using some methods similar as before) 
than the integral itself. On the other hand, once we have obtained a result, we can compute 
its average, so the comparison leads to an interesting check of the calculation. Let us see on 
the example of Gl how one performs this angular average. From ()9.5a)l we write 

in which we commuted the angular average (where dQ,\ 2 denotes the solid angle element 
in the direction n 12 ) with the integral sign and the terms depending only on y x . This is 
correct because yi is kept fixed in the process; for instance, the average of y 2 is yi, which 
is obtained by writing y 2 = y 1 — r 12 n 12 and averaging over n 12 with fixed r 12 and yi. In 
practice, computing the average (J9.29)) is not too complicated because the average of a g is 
rather simple, 
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[A more complicated example, that was useful for us, is 
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According to f)9.30|) . we must split the integration over <i 3 x into two "near-zone" and "far- 
zone" contributions, 
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The finite part at B = is necessary only for the far- zone integral. Both integrals in (19.32)1 
are now evaluated using standard methods. In the case I = 2 we find 

dni2 G il =y^ > (\nn 2 + ±-) . (9.33) 
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This is in agreement with the average of G^ computed directly with the result calculated 
from ()9.13|) or ([9.14)1 . This method of averaging has been applied for checking many other 
integrals. Even, in several cases, the method has been employed in order to determine 
some unknown coefficients. However, for this purpose the method is less powerful that 
the method of contraction with the vector yi 2 , since the latter method yields in general a 
redundant determination of the coefficients. 

The complete list of the results for the elementary integrals is as follows. 
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Inserting these elementary integrals into the expressions of non-compact terms [see (J9.3|) ]. 
and reducing to the case of circular orbits, we obtain the results reported in Appendix A. 



X. POINT-MASS REGULARIZATION AMBIGUITIES 

The computation of the multipole moments we performed so far has been carried out 
with standard techniques: standard Hadamard regularization [see Section V], and Schwartz 
distributions [see e.g. Eqs. ()9.2|) ]. The result we obtained depends on three arbitrary 
constants: the two Hadamard regularization constants u\ and u 2 introduced in Eq. ()5.3|) . 
and the constant r entering the definition of the source multipole moments through the 
analytic-continuation factor |x| B = |x/Yo| B [see Eqs. (|2.5|) ]. The constant is not a 
problem since we know that in this formalism the multipole expansion of the field exterior 



to any source is actually independent of r 2l(. Indeed we shall check in Section XII that 



ro disappears from the final expression of the energy flux (the constant ro in the source 
moments is cancelled by the same constant present in the contribution of "tails of tails" 
in the wave zone; see Eq. (jll.8j) below). However, it will turn out that the constants u\ 
and «2, which encode some arbitrariness of the Hadamard regularization, lead a priori to 
two undetermined purely numerical parameters in the expression of the 3PN quadrupole 
moment. In addition, we shall argue that because of some delicate problems linked with 
the use of the Hadamard regularization at the 3PN order, we should consider a priori a 
third undetermined parameter in the quadrupole moment. However, the important point 
is that these three parameters combine to yield one and only one undetermined constant, 
that we shall call 9, in the third time-derivative of the moment which is needed to compute 
the physical energy flux for circular orbits. Furthermore, we shall find that the constant 
9 enters the energy flux at the same level as the constant A coming from the equations of 
motion (see below), so that the energy flux depends in fine merely on one combination of 9 
and A. 

The equations of motion of compact objects at the 3PN order have been investigated 
using the ADM-Hamiltonian formulation of general relativity 0,13; an d by integrating the 
field equations in harmonic coordinates ji3,|4^|. In both approaches the compact objects are 
modelled by point-like particles described by delta-functions, and the self-field of the particles 
is removed by a Hadamard regularization. It was shown that the regularization permits the 
determination of the full equations of motion at the 3PN order except for one undetermined 
coefficient, A in the harmonic-coordinate approach and astatic i n the ADM-Hamiltonian. 
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Very likely the unknown coefficient accounts for a physical incompleteness of the point-mass 



regularization. Actually two unknown coefficients were originally introduced in [3J, [35j , but 
one of them was shown later j^J IHtJ to be fixed to a unique value by requiring, in an ad 
hoc manner, the global Poincare invariance of the Hamiltonian. On the other hand, in the 
harmonic-coordinate approach 

0,13 a new Hadamard-type regularization was developed in 
order to account for the mathematical ambiguities of the standard Hadamard regularization 
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42] . A characteristic of this regularization is the systematic use of a theory of generalized 
functions. The regularization is defined in a Lorentz-invariant way, but was ultimately shown 
to yield incomplete results for the equations of motion, in the sense that there remained the 
unknown numerical coefficient A. The complete physical equivalence between the harmonic- 
coordinate 0, 12] and ADM-Hamiltonian 0, 0, 0, E3| formalisms has been established 



38|,|44|. Indeed a unique "contact" transformation of the particles motion which changes the 
harmonic-coordinate Lagrangian (as given in Ref. [44]) into the ADM-Hamiltonian obtained 
in Ref. [37} exists. The equivalence holds if and only if the harmonic-coordinate constant A 
is related to the ADM-Hamiltonian static ambiguity by 

3 1987 
A - --c,; static - — . (10.1) 

Recently, the value ecstatic = has been obtained by means of a different regularization 
(dimensional) within the ADM-Hamiltonian approach j^. This result would mean that 
A = —1987/3080. Note that a feature of the harmonic-coordinate equations of motion 
derived in |40j, |43|, |44 1 is the dependence, in addition to A, on two arbitrary constants r\ and 



r' 2 parametrizing some logarithmic terms. However, contrary to A which is a true physical 
ambiguity, the constants r'\ and r'2 can be removed by a coordinate transformation and 
therefore represent merely some unphysical gauge constants. For instance these constants 
cancel out in the center-of-mass invariant energy of circular binaries [4oj | . 



A. Hadamard-regularization constants 

The first problem in the present calculation lies in the a priori unknown relation between 
the Hadamard regularization constants u\ and u<i introduced by Eqs. ()5.3|) and the two gauge 
constants r\ and r' 2 which parametrize the harmonic-coordinate equations of motion. Let 
us investigate more precisely the dependence of the quadrupole moment on the constants 
Ui and «2- Inspection of our computation shows that these constants come only from the 
cubic and non-compact terms obtained in Sections VIII and IX. More precisely, we find that 
the whole computation depends on ui, U2 only through the elementary integrals Y^~ 3 '°) and 
Sj~ 5 '°\ which parametrize the "self" parts, proportional to m\ or ml, of the cubic and non- 
compact terms (recall also that 5,0 ' ) is zero). See for instance the expressions ()9.3j) of NC 
terms. The relevant Y£ 3,0 ' ) and 5 '°' ) were obtained in ()8.15aj) and (J9.7|) . The dependence 
on -ui and U2 therein is 



(-3,0) 

Y- 



+ 



(10.2a) 
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= 1 4^(-)yf j> + --- ■ (io.2i» 



3 V r o, 

The dots indicate the terms independent of u± and u 2 . We take all the cubic and NC terms 
given by (|8.2j) and (|9.3|) [only the mass quadrupole is to be considered], plug into them 
the results ()10.2j) and find after summation the following part of the quadrupole moments 
depending on these constants (for general orbits): 

7 -[«i,^]= -— — i^lnM afyf + 1^2 +■•• . (10.3) 



3 c 6 \ r o . 

By Iij[ui, u 2 ) we mean the quadrupole obtained from summing all the terms computed in the 
previous sections, i.e. depending on the Hadamard-regularization constants u\, w 2 (as well 
as, of course, the constant tq). On the other hand, we found that many of the "interaction" 
terms, proportional to m\m2 or mim|, depend on time-dependent logarithms of the ratio 
?\2 = ru/ro, where r$ is the constant dealing with the behaviour of the moments at infinity. 
See for instance the elementary integrals (|9.8jl . The effect of the result (jlO.3)) is to "replace" 
a part of the latter logarithms of f\ 2 by some corresponding logarithms of the ratio r^/ux 
(and ditto with m 2 ). The remaining logarithms stay as they are as logarithms of the ratio 
f 12 . Thus we can re- write the dependence of the quadrupole on U\ and U2 through the 
logarithms of r^/ui and r 12 /-u 2 in the form 

I ij [m,U2] = (^^ln( r -^] a?yi> + l~2)+... . (10.4) 



C 6 \ Mi 



All the other logarithms, present in the dots of Eq. ()10.4j) . are of the type In in 
this paper we assumed nothing about the values of u\ and u 2 . In particular we did not 
assume any relation between m, w 2 and the gauge constants r\, r' 2 that parametrize the 
final equations of motion in harmonic coordinates jiol. l43^. However, when computing the 
energy flux we shall need to obtain the third time-derivative of the quadrupole moment, 
and for that purpose we shall replace the accelerations by their expressions obtained from 
the 3PN equations of motion, depending on r\, r' 2 . As a result the third time-derivative 
of the moment will depend on u\, w 2 as well as on r\, r' 2 . Therefore, we definitely need to 
control the relation between u±, w 2 and r\, r' 2 ; then we shall have the quadrupole moment 
expressed solely in terms of r\ and r' 2 and we shall check that the latter constants can be 
removed by the same coordinate transformation as in the equations of motion, and thus 
that the final expression of the physical energy flux must be independent of these constants. 
From Eq. ()10.4|) we can write 

i^-K u 2 ] = hjir'u r' 2 ] + In (11) afyt + 1^2. (10.5) 

The notation for Iij[r\, r'2] is clear: we mean the sum of all the contributions obtained in 
the previous sections, but computed with r\, r' 2 in place of the regularization constants u\, 

We shall now look for the most general In {^^j that is allowed by physical requirements. 
In this connection recall the spirit of the regularization: the constants u\ and u 2 reflect some 
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incompleteness of the process, that may or may not be fixed in a given computation, and 
therefore they should be kept completely arbitrary unless there are some physical arguments 
to restrict their form. In particular, when used in different computations, these regularization 
constants have no reason a priori to be the same. For instance, in the present computation 
of the moments, the constants u\ and w 2 are a priori different from the constants s\ and s 2 
which were originally used in the 3PN equations of motion (see Eq. (2.3) in 43]). They are 
a fortiori different from the constants r\ and r' 2 chosen to parametrize the final equations 
of motion (Eq. (7.16) in ji^|). See also the discussion in Section VII in Ref. where 
we determined the general form of the relation between s±, s 2 and r\, r'2 by imposing the 
polynomial mass dependence of the equations of motion, the correct perturbative limit, and 
the existence of a conserved energy. Here we shall basically do the same in order to restrict 
the form of the relation between U2 and r\, r' 2 . Note that a priori the logarithms 
ln(r'i/wi) and ln(r' 2 /w 2 ) can depend on the masses mi and m 2 . To determine just what 
combination of masses is allowed we make (similarly to the equations of motion) two physical 
requirements: (i) that the quadrupole moment be a polynomial function of the two masses 
mi, ^2 when taken separately; (ii) that the perturbative limit (corresponding to v — > 0) 
not be affected by this possible dependence over the masses. Because of the factor m\ in 
front of the log-term in ()10.5|) . and because the acceleration a\ brings another factor m 2 , 
the most general solution for this logarithm in order to satisfy the requirement (i) is to be 
composed of: a pure numerical constant (say £), plus a pure constant (say k) times the mass 
ratio m/mi, plus a constant times m/m 2 , next five terms involving the mass ratios m 2 /m 2 , 
m 2 /mi/m 2 , m 3 /mf, m 3 /m\/m2 and m 4 /m 3 /m 2 . Each of these terms must be such that it 
does not violate the perturbative limit [our requirement (ii)]. This means that they should 
involve, in a center-of-mass frame, a factor v 2 at least. We readily find that the only two 
admissible terms in this respect are the first two in the previous list (with constants £ and 
k). So we end up with the most general admissible solution 



where £ and k denote some arbitrary purely numerical constants (for instance rational 
fractions). This result is similar to the one obtained in the 3PN equations of motion, 
concerning the relation between si, s 2 and r\, r' 2 . See Eqs. (7.9) in Ref. [43], where the 
determination of the constant analogous to £ was possible from the requirement of existence 
of a conserved energy (and Lagrangian) for the equations of motion. 

We now check that the logarithms of ri 2 /r\ and r 12 /r' 2 in the quadrupole moment, which 
are of the form 



can be eliminated by the same coordinate transformation as found in Ref. 43] for the 
logarithms in the harmonic-coordinate equations of motion. [As concerns the logarithms of 
r i2/ r o i n the moment they cannot be eliminated by a change of coordinates but will match 
precisely with corresponding logarithms present in the "tails of tails" at infinity] We look 




(and idem with 1 <-> 2) , 



(10.6) 




(10.7) 
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for a coordinate change of the type considered in Section VI. A of [43]: namely Sx^ = 
where £ M = 77^ is a 3PN gauge vector given by 

f „ = ^! a „(fl + li) . (10 . 8) 

We have factorized out m 3 (where m = mi + m 2 ) so that t\ and e 2 , which are constants 
or mere functions of time t, be dimensionless. The corresponding change of the particle's 
trajectories is given to this order by the regularized value of the gauge vector at the location 
of the particle (see Section VI. A in We obtain 

G 3 7Tl 3 

kv\ = - £ 2-^y l i2 , (io.9a) 
c '12 

%2 = ei-^TVL. • (10.9b) 



12 



Since the quadrupole moment starts at the Newtonian level with the usual m\y^ %3> + 1 
we easily find its coordinate change as 



S^Iy = 2m 1 ?/ 1 < %/ 1 > + 1^2 
G 3 m 3 
c 6 r 



-2m 1 e 2 -^y< i y^ + 1^2. (10.10) 
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By comparing this with Eq. ()10.7j) (using the Newtonian particles acceleration), we find 
that the gauge transformation required to eliminate the logarithms is 



22 m?m 2 , ( T\i \ 
e 2 = --^/lnU , 10.11a 



£ i = ---^ ln ( — ) ' ( 10 - llb ) 



m° \ r 



in complete agreement with Eq (7.2) in Ref. |43j. In summary, not only these logarithms 
disappear when considering physical quantities associated with the equations of motion (such 
as the invariant energy), but they will also cancel from physical quantities associated with 
the wave field at infinity, viz the invariant energy flux we compute in Section XII. 



B. Special features of the regularization 

We now discuss some subtleties of the Hadamard regularization which motivate the in- 
troduction in the quadrupole moment, in addition to £ and k considered in Eq. ()10.6|) . of 
still another constant (however, see below for the definition of a single constant 9) . 

N on- distributivity of the Hadamard partie finie. By "non-distributivity" we mean the fact 
that the regularization of a product of two functions F and G, singular in the sense of (|5.1jh 
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does not equal, in general, the product of the regularized functions: (FG)i ^ (F)i{G)\. For 
instance, with U = Gm\jr\ + Gm^jr^ the Newtonian potential, we have {U n )\ = [(U)i] n 
for n = 1,2,3, but (£/ 4 )i = [(t^)i] 4 + 2[(£7)i] 2 [([/) 2 ] 2 . An immediate consequence is that 
the product of a singular function F with a delta-function does not equal, in general, the 
product of its regularized value with the delta-function: F8\ ^ {F)i8\. Here we are assuming 
that the three-dimensional integral of the product of F with 8\ = 5(x — yi) gives back the 
regularized value (F)i. Notice that only at the 3PN order does the non-distributivity play 
a role. Up to the 2PN order, the distributivity holds for all the functions encountered in the 
problem (hence the computation of the moments as was done in 0| is correct). 

The non-distributivity at 3PN has an important bearing on the choice of the stress-energy 
tensor for describing point-particles. In this paper, we adopted the most naive choice for 
the stress-energy tensor. See Eq. (|5.4|) above, which is equivalent, at 3PN order, to 



= - miV ' V l g(x - y x ) + 1 <-> 2 . (10.12) 



for 

y/{99pa)iv{vl/(? 

Namely, we assumed that the whole factor of the delta-function consists of a regularized 
value at point 1. But because F8\ ^ (F)i8\ 1 we could obtain a different result by choosing 
another stress-energy tensor, defined by replacing the factor of the delta-function in (|10.12|h 
or part of it, by a function depending on any field point x and such that its regularized 
value when x — > y x is the same. In fact, a specific form of the stress-energy tensor of 
point-particles, compatible with the Hadamard regularization, was advocated in Ref. ji^ 
and used to compute the 3PN equations of motion This form, given by Eq. (5.11) in 
Ref. H, reads 



T- = miV ^ 4*ZgL + 1^2. (10.13) 



Choosing one or the other form of stress-energy tensor does make a difference in our com- 
putation. Consider for instance the term SI(1) = J d z ~x.x l ^a. We find that the result for 
this term, when computed using the tensor 1)10.13)1 i.e. using a = T 00 + T M , differs from the 
original result by the amount 



Q2 3 



2 1 

Z n <\p> _ ± 7 ,< i 7 ,J > 
2«i 2/i g u i u i 



+ 1^2. (10.14) 



There is also a modification ASII(l) but which is of the same structure (with different 
numerical coefficients). 

On the other hand, some terms in our computation would be different if the regularization 
would be distributive. For instance, if for computing the term SI(16NC) we take into account 
the non-distributivity (as we did), we find the result ()6.9|) . namely 

SI(16NC) = (x <l x j> U^ c) ^ + 1^2. (10.15) 

If instead we incorrectly assume that the partie finie is distributive, then we get 

Sl(i6NC) distr = ^yfyfvf (uS c) ) i + 1 - 2 . (10.16) 
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The difference between the two results is not zero: 



ASI(16NC) 



2 G 2 m\ 
15 c 6 



vfvf + 1 



(10.17) 



The same happens with the other terms VI(IONC) and VI(12NC); each time the structure 
of the difference is the same as in ()10.14|) or ()10.17j) . 

Violation of the Leibniz rule by the distributional derivative. In Ref. 



41 a new kind of 



distributional derivative of singular functions of the type F was introduced. It was found 
that it is impossible to define a derivative satisfying the Leibniz rule for the derivation of the 
product, but that a mathematical structure exists when we replace the Leibniz rule by the 
weaker rule of "integration by parts". The latter rule can be seen as an integrated version 
of the Leibniz rule (see Section VILA in 
derivatives were proposed in 
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More precisely, two different distributional 
a "particular" derivative, and a "correct" one. Both 
derivatives reduce to the derivative of the standard distribution theory when applied to 
smooth test functions with compact support. The particular derivative is simpler to use in 
practical computations, but the correct one is more satisfying because successive derivatives 
to any order commute. 

Previously we performed numerous simplifications, with the help of the Leibniz rule, to 
arrive at the form of multipole moments given by Eqs. ([4.2)1 . Thus we made some errors 
because of the violation of the Leibniz rule by the distributional derivative. The strategy 
adopted in Ref. 43J] was to keep track of all these error terms and to compute them using the 
particular and correct derivatives of In the present paper we shall proceed differently. 
We simply give an example. When simplifying the moment to arrive at the simple-looking 
term SI(39) in Eq. (|4.2aj) . we "forgot" to include the error term 



ASI(39) = FP [ rf 3 x|x| B ^- 7 - \A(U 4 ) - AU 3 AU - 12U 2 d a Ud a U] . (10.18) 

3nGc b b=o J 

Clearly this term would be zero for any derivative satisfying the Leibniz rule (in a distribu- 
tional sense). However, computing it by means of for instance the "particular" derivative 
(defined by Eq. (7.7) in Ref. 



we find that it is not zero, but 



64 G 2 m\ 



ASI(39) = + (10.19) 



Again this result has the same type of structure as found previously. We have checked that 
all the terms coming from the violation of the Leibniz rule have the same structure, either 
of type m\af l y 3 ] > like in ()10.19|) or of type m\vf l v { > . 

Cubically non-linear self-interaction terms. We take the example of the self contribution in 
the term SI(35NC). This term is 

SI(35NC) = ^— FP [ d 3 x\±\ B x L zlf C) dlU . (10.20) 

irGc 6 b=o J 3 3 ' 
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The "self" part of this term corresponds to that part in zf^ C ^ which is proportional to m\, 
in the sense that 



HNC) 



G' z m'i <j at dj) In r x + ~-v{ dij In n + — 5^< m In n 

1 1 OZ 1 



km 



32 r? 



(10.21) 



[see Eq. (I9.1d|) ]. and that part of U due to 1 itself, i.e. U = Gm\jr\ + C(m 2 ). The resulting 
term, proportional to m\, is ill-defined in distribution theory because the delta-function, 
coming from the distributional derivative of l/r\ as given by Eq. (|9.2b[) . is multiplied by the 
terms in (jl0.21|) which are singular at point 1. The partie finie pseudo-functions and their 
derivatives proposed in Ref. |4l| permit us to give a mathematical meaning to such ill-defined 
terms. The "particular" derivative of l/r 1 reads d^{l/rx) = <9| ; (l/r 1 ) ordinary + D^l/ri], 
where the purely distributional part is 



D 



1 



Att ( „■„■ 15 ^aa 
3 V 2 1 



Si ■ 



(10.22) 



[Compare this with the result (J9.2bj) of distribution theory] We easily compute the effect of 
this new derivative on the self part of the term (jl0.20j) . Once again we find the same type 
of structure as before: 



ASI(35NC) 



G 2 m\ 



64 <t 3> 38 
3 i fi 15 1 1 



u n < V > 



10.23) 



Similarly we checked that all other self-interaction contributions take the same form with 
simply different numerical coefficients. 



C. Definition of the ^-ambiguity 

As we have seen the structure of the possible terms associated with the previous subtleties 
in the Hadamard regularization is limited to only two types, either mfa^y^ or mfffV^. 
The first type was already considered in Eqs. (|10.5p - ()10.6jl . where it yielded the arbitrary 
constant £. Thus, modulo a redefinition of £, we do not need to consider this term. The 
other type, given by m^v^v^, was not considered earlier. Therefore, motivated by the 
previous discussion, we shall from now on add such a term to the multipole moment, with 
a new constant in front, say (. In summary, we consider three types of "ambiguous" terms 
(in the sense of ji^, EB|), parametrized by the two constants £, k of Eq. ()10.5|) and the (. 
The quadrupole moment we finally consider in this paper is thus 

I lJ =I ij [r' 1 ,r' 2 ] + AI lJ , (10.24) 

where ijjfr'i, rV| denotes the computation we have done in Sections VI-IX (i.e. the sum 
of all the terms, defined for general orbits, and given for circular orbits in Appendix A), 
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when expressed by means of the same regularization constants r\, r'2 as the ones appearing 
in the 3PN equations of motion (we know that these constants are pure gauge). Now the 
undetermined part reads as 



44 G 2 m 3 
In a center-of-mass frame we get 



mi J 



+ 1^2. (10.25) 



44 G 2 m?v 2 r 

AI tj = ? ~ 6 [(£ + 2k) a<V> + C u < V > ] (10.26) 
o c 

(where x l = y\ — y l 2 , v % = dx % jdt and a* = dv l /dt). The constants £, /t and £ will be 
left unspecified in the present paper. It could be possible that the more sophisticated 



regularization procedure of |41l . I42j determines some of these constants. However, the point 
for our purpose is that we are going to show that the physical energy flux for circular orbits 
depends only on one parameter. Indeed, the flux depends on the third time-derivative of 
the quadrupole, and by a straightforward computation (using the Newtonian equations of 
motion) we find that, in the case of circular orbits, the third time-derivative of ()10.26j) is 



Af = ^^7 3 ^<V>, (10.27) 



352 Gm 2 v 2 
3 r 3 

where # = £ + 2fi;-|-£isa single unknown constant. Therefore, the ambiguous part of the 
physical 3PN flux, as concerns this effect, depends in fact only on 9. It is given (for circular 
orbits) by 

^ = i/,?'A^ = ^V{-f^}. (10.28) 

In addition to 9, the flux will depend also on the constant A coming from the equations 
of motion 0,0. However, we shall find that, in the case of circular orbits, both 9 and 
A enter the flux at the same level, so the flux depends only on one combination of these 
constants: A — |#, from the end result ()12.9j) below. Further work, supplementing the 
Hadamard self-field regularization by suitable extensions and alternative methods, may be 
required to determine the constants 9 and A. 

XI. THE BINARY'S MULTIPOLE MOMENTS 

The computation of the moments is now almost complete. The remaining terms are 

(i) the "odd" terms: SI(ll), SI(12). These terms involve the fifth (odd) power of 1/c (2.5PN 
order). They appear because of the expansion of retardations in the potentials ()3.6|) : they 
are pure functions of time, parametrized by Qij(t) and Q(t) [see (13. 9|) ]. The sum of the 
two odd terms has been computed in the equations (4.9) and (4.12) of Ref. 
present notation, in the quadrupole case, we have 



26J. With the 



si(ii) + si(i2) = ^ |-|gg 4 Q i>fc - yQ%>Q ~ ^Q<i3>Q m } ■ (n-i) 
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These terms do not contribute to the flux for circular orbits. 



(ii) the "divergence" terms: SI(2), SI(8), SI(9), SI(10), SI(39), SI(40), SI(41), SI(42), SI(43), 
SI(44), SI(45), SI(46), SI(47), SI(48), SI(49), SI(50), SII(8), SII(IO), SII(12), SII(14), SIII(3), 
VI(6), VI(30), VI(31), VI(32), VI(33), VI(34), VI(35), VII(7), Tl(9). The integrand of these 
terms is made of the product of |x| B and a pure divergence diA or AA. Their computa- 
tion makes use of the same techniques as those employed in Sections VII-IX, but with the 
notable simplification that because of the divergence one can perform an integration by 
parts, and that as a result the elementary integrals contain explicitly a factor B (due to the 
differentiation of |x| s ) so their computation is quite easy. See the results in Appendix A. 

(iii) four particular terms that we have left out because their sum is in fact zero: 

SI(22) + SI(23) + SI(32) + SII(ll) = . (11.2) 

We sum up all the terms given in Appendix A, plus the undetermined correction given by 
f)l().2f)jl . and obtain the expressions of the 3PN mass-type quadrupole moment, 2PN mass- 
type octupole moment and 2PN current-type quadrupole moment of the compact binary 
moving on a circular orbit. (Note that most of the investigation of this paper is valid for 
general orbits, but we are interested in inspiralling binaries whose orbit is quickly circularized 
by radiation reaction.) The 3PN mass quadrupole reads 

(r 2 48 t \ 

Axij + B —Vij + — -x^Vfrvf J + 0(7) , (11.3) 

where the third term is the 2.5PN odd term, and where 
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The mass parameters are: m = mi + 777-2, $ m — m i — m 2, A 4 = ^1777.2/777 and v = /x/m. The 
post-Newtonian parameter is 7 = Gm/(rc 2 ) = 0(2) [see Eq. (J5.8j) ]. The logarithms depend 
either on the constant ro associated with the finite part at infinity (recall |x| B = |x/ro| B ) or 
on the "logarithmic barycenter" r'o of the regularization constants r\ and r'2 (see Section 
X), defined by 777 In r'o = mi lnr'i +m2 In r'2. We shall investigate in Section XII the fate of 
these constants r and r' . In addition the moment depends on the unknown constants £, K 
and ( introduced in Eq. ()10.26|) . The 2PN mass-octupole and 2PN current-quadrupole are 
free of any of such constants and given by 
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The higher multipole moments which are needed in the 3PN energy flux are the 1PN current 
octupole, 1PN mass 2 4 -pole, Newtonian current 2 4 -pole and Newtonian mass 2 5 -pole. For 
these moments we simply report the expressions already obtained in Ref. jl^ . 
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As proved in Refs. |21l . I26j the multipole moments II and Jl are not the only source 
moments entering the radiation field. However, the other moments, denoted Wl, Xl, Yl 
and Zl, parametrize a (linearized) gauge transformation in the exterior field, and as a result 
make a contribution to the non-linear radiation field at a quite high post-Newtonian order: 
2.5PN. It is always possible to re-express the radiation field in terms of solely two sets of 
moments, denoted Ml and Sl, given by some non-linear functionals of the moments II, Jl, 
Wl, Xj„ Yt. and Z L , but differing from I L and Jl starting at the 2J3PN order (see Section 
VI in 

Mij is related to by 



for a discussion). From the equations (4.20)-(4.24) in [26j] the 3PN quadrupole 
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For the other moments there is no correction to be made at this order [for instance Sij = 
Jij + 0(5)]. Actually we observe that W is zero for circular orbits, and thus we shall from 
now on replace all the moments II and Jl by the corresponding M L and Sl- 
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Finally we need to relate the moments Ml, Sl to the "radiative" moments, say Ul (mass- 
type) and Vl (current-type), which play the role of observables associated with the radiation 
field at infinity. Since such a relation has already been worked out at the 3.5PN level in 
Ref. 27[, we simply report the main result, which concerns the mass-quadrupole radiative 
moment £7^, that is 
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This formula is valid through 3.5PN order, modulo the odd-order 2.5PN and 3.5PN terms 
that we do not show because they do not contribute to the flux for circular orbits. The only 
contributions coming from Eq. ()11.8|) are the 1.5PN tail, and 3PN "tail of tail" integrals. 
In the flux we shall derive below the terms at the orders 2.5PN and 3.5PN are due to the 
tail integrals in higher multipole moments (see Ref. 0] for details). 



XII. THE ENERGY FLUX OF CIRCULAR COMPACT BINARIES 



For general sources, the total energy flux (or gravitational luminosity C) to the 3PN 
order is composed of an "instantaneous" contribution — i.e. a functional of the multipole 
moments Mr an d Sl at the same instant — , and a "tail" contribution. We shall now follow 
the study in 27[ of the occurence of non-linear effects in £ up to 3.5PN order. Following the 
equation (4.18) in Ref. |2j| we split C into an instantaneous part, a tail part, a tail square 
part, and a tail of tail part: 



£ — Anst + Aail + £(tail) 2 + Aail(tail) • (12.1) 

As all the parts involving tails have already been computed for circular binaries j27j . we 
need only to compute the instantaneous part which is given by 
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(see e.g. (4.15) in [27]). To obtain £ inst we compute the time derivatives of the multipole 
moments. At this stage we need a new input, namely the 3PN equations of motion of circular 
binaries which are crucial in the differentiation of the 3PN quadrupole moment. As recently 
obtained 0,0 the 3PN orbital frequency of the circular motion reads as 
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The inverse of this formula gives the post-Newtonian parameter 7 as a function of the 
frequency related parameter x = (Gmco / c 3 ) 2 / 3 , 
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Note that (|T2~3|) or (jTOjl involve the same constant r'o as in the 3PN mass quadrupole 
moment (fTOjl - ffTOjl . 

Taking all the expressions of the multipole moments found in Section XI, computing 
their time- derivatives according to the latter circular-orbit 3PN equations of motion, and 
inserting them into ()12.2|) we then arrive at the following instantaneous part of the flux, 
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where we recall that = £ + 2k + £. Next we simply add the known other contributions. 
The tail one is due to such terms as the 1.5PN integral appearing in Eq. (J11.8)) [and other 
equations corresponding to hig her multipole moments]. The result is derived to the 3.5PN 
order in Eq. (5.5a) in Ref. [271 ]: 
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Second, the tail of tail comes from the 3PN term in Eq. (jll.8j) . and the tail square from the 
square of the 1.5PN term. The sum of these parts reads, following Eq. (5.9) in |27| . 
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where C denotes the Euler constant (C = 0.577 • • • ), and where the constant ro is the same 
as the tq occuring in the mass quadrupole moment (jll.3j) - ()11.4j) . Thus, the energy flux, 
complete up to the 3.5PN order, reads 
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We observe that the constants r have cancelled out between the instantaneous flux £ inst 
and the part £(taii) 2 +taii(taii)- This cancellation is to be expected for any source: see a proof 
(Eqs. (4.14) there) where it is shown that the tails of tails at the 3PN order depend 

^J — v ' To5 ' u -Jr- M^f . Using our 
(9(5) we find that indeed the r 's cancel out. The fact 
that we have recovered the expected dependence on r of the source quadrupole moment is 
a good check of the computation. 

On the other hand, the point-mass regularization constant r' still remains in the flux 
()12.8j) . This is because the energy flux is not yet expressed in a coordinate- independent 
way, as the post-Newtonian parameter 7 depends on the distance between the masses in 
harmonic coordinates. To find a truly coordinate-independent result we must replace 7 by 
its expression given by (|12.4j) in terms of the frequency-related parameter x. With this 
change of variable, at long last we obtain our end result: 
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In the above expression the constant r'o has cleanly disappeared. Of course, this was to 
be expected because we have seen that r' is pure-gauge; nevertheless this cancellation 
constitutes a satisfactory test of the algebra. However, the result still depends on one 
physical undetermined numerical coefficient, which is a linear combination of the equation- 
of-motion-related constant A and the multipole-moment-related constant 9. On the other 
hand, our final expression ()12.9|) is in perfect agreement, in the test-mass limit v — ► 0, 
with the result of black-hole perturbation theory which is already known to a very high 
post-Newtonian order 46[ 47]. 
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APPENDIX A: RESULTS FOR ALL THE TERMS 



For the mass quadrupole we factorize out a factor /x = mv in front of all the terms. We 
denote v % = yj ^w 1 , so for instance Wij = jj^v^v^ (and ctij = x^x^). For the current 
quadrupole all the terms have to be multiplied by 5m/m L < iXj > , where 5m = m\ — rri2 
and Li = fi EijkX^v k is the angular momentum. For the mass octupole we factorize out 
n5m/m = u5m. For simplicity the constants tq and r'o in the logarithms are set to one. 
In the case of the 3PN mass quadrupole, to the sum of all these terms one must add the 
undetermined contribution given by Eq. f)10.26|) in the text. 



1. The 3PN mass quadrupole 

Miscellaneous: 



SI(22 + 23 + 32) + SII(ll) = , 

VI(16 + 20) + VII(6) = VI(19) 



Compact term at Newtonian order: 



SI(1) 



Compact terms at 1PN: 



1 + -(1 



5v) 



7 



r 



+ — (149- 573z/ + 354z/ 2 
16 V 



13-61z/ + 25z/ : 
29u 3 ) x 



(Ala) 
(Alb) 



(A2a) 



SII(l) = — 7 [(-8 + 24z/ + 7(20 - 52z/-20z/ 2 ^ 
56 

+ 7 2 (-23 - \lv + 160z/ 2 - 55u 3 ))xij 

+ (8 - 24z/ + 7 (4 - 28z/ + 44z/ 2 ) 

+ 7 2 (-13-9z/ + 238z/ 2 + 35z/ 3 ))i%] , 

VI(1) = ^7 [(-8 + 24z/ + 7 (28 - 92z/ + 20z/ 2 ) 

+ 7 2 (-75 + 259z/- 176z/ 2 + 13z/ 3 ))% 

+ (8 - 24z/ + 7 (-4 + Ylv + 4z/ 2 ) 

+ 7 2 (15 - 157z/ + 414z/ 2 + 7z/ 3 ))w ii ] . 



(A3a) 



(A3b) 



57 



Compact terms at 2PN: 
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VII(4) = T J^7 3 (3 - 2u - 34z/ 2 + 8v 3 )(w tJ - x l3 ) , (A8f) 

VII(5) = "t|7 3 (-1 + 61/ - 18z/ 2 + 8z/ 3 )K- - , (A8g) 

111(2) = ^7 3 (8 - 67z/ + 134z/ 2 - 12z/ 3 )(^, - %) . (A8h) 



T-terms at 3PN: 



SI(17) = -U 3 z/ [(7-Wv)wij+ (-15 + 36z/ + 54z/ 2 )%] , (A9a) 
9 
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SI(18) 
SI(19C) 

SI(21C) 
SI(24) 
SI(25) 
SII(9) 
VI(14) 
VI(15) 
VI(17) 
VI(18) 
TI(5) 



Cubic terms: 



- 7 3 z/ [-3^ + (-1 - I2u + \Su 2 )x i3 ] , 



= r 



12 
1 



18 



--(1 - v)wij + —(2 -7u+ 70u 2 + 12z/ 3 ) 

7 3 z/ [(-2 - %v)wij + (-1 + 22v + \2u 2 )x i3 ] , 
7 3 z/ [llwij + (-15 + 36z/ + 54z/ 2 )%] , 



9 



■- 7 3 z/ [-i&tf + (-3 - 36z/ + 54z/ 2 )%] , 
7 3 (l-9^)(-^ + %), 
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7 3 (20 - l\v + 3Qu 3 )(w i3 - x^) , 



-— 7 3 z/(17 - 98z/ + 36z/ 2 )(u^- - x l3 ) , 
^7 3 z/(-43 + 94z/ + 36v 2 )(w i3 - , 



42 
1 



7 3 (4 - 3v - 72v 2 + 36z/ 3 )(t% - Xij) , 



189 



7 3 (4 - 9u - 2\v 2 + 13u 3 )(w i3 - x^) 



SI(26 + 27 + 28 + 29 + 30 + 34 + 36) 



7 3 z/ 



74 152 

v I ajji + ■ 

15 3 ) 13 V 15 



128 , 1024 \ 
mr w 



225 7 



SII(13) 



16 
225 



7 3 ((45 - 75i/) lnr - 9 - 25z/)(-% + w tj ) , 



VI(22 + 23 + 24) = — 7 3 (601nr- 7- 265z/ + 275i/ 2 )(£ ii -w ii ) . 

525 



Non-compacts terms: 



SI(5NC) = --7 2 (4+10z/ + 7 (2-6z/-47z/ 2 ))x lj , 



SI(19NC) = 
SI(20) = 
+ 

SI(21NC) = 



1 71 



7 



6 36 



1 77 



77 + - n V ) 7 + T^^lnr — - - —=-V 7 Wij , 



1 1 15 
77 25 



6 450 



4 4 . 

1 — v mr 

15 3 ) 225 6 



-^-v + 6z/ 2 ) 7 3 % 



4 4 \ , 77 77 , , 

1 v mr H v 7 Wi, , 

15 15 / 225 225 1 1 3 ' 



8 



26 



1 



— — — —v ) lnr + — — — — v — v 2 I 7 3 a;i 1 
15 3 225 36 ' ' 3 
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f f 8 8 A , 26 397 \ ,„ , Ai1 

+ 1 v lnr 1 v Ywu , (Alld) 

VV 15 15 J 225 450 J 1 3 ' V ; 

SI(33NC) = (-^+f- 2 )7%+(-^lnr-ig^ 7 3 ^, (Alle) 

nx.or.T^x /Y 16 64 V 1° 28 II 48 I 45 2 \ 3~ 
SI(35NC) = ((-- - -u) In, + — + _„ - 

( (is 16 A, 1028 868 "\ 3~ /.., r > 

+ ((i5-yT nr -i2r + ^T' w '" (Allf) 

SI(37NC) = (~v - "v 2 ) + (-^flnr + 7 3 »a , (Allg) 

SI(38NC) = {^v + |^ T % + (^Inr + 7 'u, y , (Allh) 

SII ' 4NC » ^ ((I - r) lnr + 1S5 - 315 " + f" 2 ) ^ - *«) • < An " 

, fT/nrlin , // 48 8 \ 968 10792 236 ,\ 3 _ . . /A11 .. 
VI(25NC) = ((-- + -„J lnr- — + — „ - T 3 (*« - Wy ) , (A1U) 

VI(26NC) + VI(27NC) + VI(28NC) + VI(29NC) = ( ( ^- - ^-is) lnr 

\\105 15 / 

2644 488 124 2 \ 3 . A A . , A111 . 

H v v l Yixij - wa) . (Allk) 

1575 175 63 j nj " K J 

Terms at 2.5PN: 



48 

SI(11 + 12) = —"j b/2 vx <i w j> . (A12) 



Divergence terms: 



SI(2) 


= , 


(Al3a) 


SI(8) 


= , 


(A13b) 


SI(9) 


= 0, 


(Al3c) 


SI(10) 


= , 


(A13d) 


SI(39) 


= 0, 


(Al3e) 


SI(40) 


= 0, 


(A13f) 


SI(41) 


= -g7 3 (-% +Wij) , 


(A13g) 


SI(42) 


= 0, 


(A13h) 


SI(43) 


= 0, 


(A13i) 


SI(44) 


1 3/ - ^ 
= g7 {-Xij + Wij) , 


(Al3j) 


SI(45) 


= 7 3 (-Xij + Wij) , 


(Al3k) 


SI(46) 


= g7 {Xij - Wij) , 


(A131) 
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SI(47) 



= 0, 



SI(48) = --^(-Xij+Wij) , 

SI(49) = 16 7 3 (% - Wij) , 

SI(50) = p^j-Wij), 

SII(8) = , 

SII(IO) = , 

SII(12) = ^ 7 3(£. ._«,..) , 
72 

SII(14) = _ 7 3(_^. + ^.), 

SIII(3) = , 

VI(6) = , 

VI(30) = , 

VI(31) = , 

VI(32) = , 

VI(33) = ^(xij-Wij), 

VI(34) = , 

VI(35) = -1 7 3(^._^..), 

VII(7) = , 

TI(9) = . 



2. The 2PN current quadrupole 



Compact term at Newtonian order: 



VI(1) = - [-8 + 4 7 + 7 2 (-15 + 88u + 3z/ 2 )] 



Compact terms at IPN: 



VI(2) = 7 [-2 + 2u + 7 (1 + v - 4z/ 2 )] , 
VI(3) = 7 z/[-2+ 7 (-l + 4z/)] , 

VII(l) = ^ [2 - Av + 7 (-7 + 16i/ - 3z/ 2 )] 
28 

TI(1) = [2 - Av + 7 (-7 + 16i/ - 3z/ 2 )l 
56 



Compact terms at 2PN: 
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VI(7) 


= 7Ml " v) , 


VI(8) 


= 2 7 2 (-l + 2z/) , 


VI(9) 


= -f 2 v(l + v) , 


VI(IOC) 


= -4 7 2 z/ , 


VlflONOl 

V 11 v. f 


3 2 
2 7 


VI(ll) 


= 2 7 V(1 - u) , 


VI(12C) 


= o, 


VI(12NC) 


= y(i-H, 


V 11 1 I 




VII(2) 


= y(l-3z/+z/ 2 ) , 


VII(3) 


= y"(l-"), 


VIII(l) 


= ^(1-3i/)(-1 + i/), 


TI(3) 




TI(4) 


= y"(l-"), 


TII(l) 


= ^(l-3,)(-l + ,). 



F-terms at 1PN: 

VI(4) = 2[-2 + 4z/ + 7 (l-8z/ 2 )] , 

VI(5) = -! 7 „[2 + 7 (l-6i/)] , 
TI(2) = . 

F-terms at 2PN: 

'Y 2 

VI(16 + 20) + VII(6) = -HI -2u- 3z/ 2 ) 

VI(19) = ^-(l-2u-3u 2 ), 

VI(21) = -^ 7 2 zy(5-6z/) , 
VI(25C) = 2 7 2 (-l + 2z/) , 
VI(26C) = ^u(l + 3z/) , 
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S-terms at 2PN: 



T-terms at 2PN: 



Cubic terms: 



Non-compact terms: 



VI(27C) 
VI(28C) 
VI(29C) 

TI(6) 
TI(7) 

TI(8) 



;i-2^), 



11 
3 

0, 

2 2 
-7 v , 

3 1 ' 



7 
2 
0, 



2_ ( _ 2 + 4z/-3z/ 2 ), 



-^(7-3,). 



VII(4) = l 7 2 (-2-3z/)(l-2z/), 



VII(5) 
TII(2) 



-^7M2 + 3z/), 
. 



:7 2 (1 +4z/-9z/ 2 ) , 



1 



VI(14) 

VI(15) = ^/(l -9//) 
VI(17) 
VI(18) 
11(5) 



i 7 2 z/(4 + 9u) , 
| 7 V(2 + 3i/) , 



1 



168 



7 2 i/(2 + 3i/) 



VI(22 + 23 + 24) = 2 7 2 



1 9 

3-r 



VI(25NC) 



7 2 l-: 



VI(26NC + 27NC + 28NC + 29NC) = 7 2 



Divergence terms: 



VI(6) = , 
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VI(30) = , (A23b) 

VI(31) = , (A23c) 

VI(32) = , (A23d) 

VI(33) = , (A23e) 

VI(34) = , (A23f) 

VI(35) = , (A23g) 

VII(7) = , (A23h) 

TI(9) = . (A23i) 



3. The 2PN mass octupole 



Compact term at Newtonian order: 



SI(1) = - [-8 + 7("4 + 16i/) + 7 2 (13 + 24z/ + 15z/ 2 )] x ijk . (A24) 

o 

Compact terms at 1PN: 

SII(l) = ^ 7 [(30 - 60z/ + 7 (-75 + 120z/ + 45z/ 2 ))x ijfc 

+ 3(-20 + 40z/ + 7 (-10 + 60z/-70z/ 2 ))a; < iW ife> ] , (A25a) 
VI(1) = — 7 [(30 - 60z/ + 7(-105 + 240z/-45z/ 2 ))x ijfc 

+ 3(-20 + 40z/ + 7 (10-20z/- 10i/ 2 ))x<^ fc> ] . (A25b) 
Compact terms at 2PN: 



SI(3) = -4 7 2 (1 - 3u + v 2 )x ijk , (A26a) 

SIII(l) = ^7 2 (1 - 3i/)(l - u)(-105x ijk + 300x <iWjk> ) , (A26b) 

VI(2) = ^7 2 (1 - 3i/ + v 2 ){lhx l]k - 30x <iWij> ) , (A26c) 

VI(3) = |:7M-1 + u)(-15x ijk + 30x <iWjk> ) , (A26d) 

VII(l) = ^7 2 (1 - 3i/)(l - u)(-105x ijk + 300x <iWjk> ) , (A26e) 

TI(1) = ^7 2 (1 - 3i/)(l - i/)(-33x ijfc + Ux <iWjk> ) . (A26f) 

F-terms at 2PN: 

SI(4) = -—-f 2 iy[(15-30iy)x ijk + 60x <iWjk> } , (A27a) 
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Oil OKj J — 


— '-v 2 \— (A ^ — (SD;y — (\V\u 2 \^r ■■, -1- "\( — 1 D -I- 4D;/W -7/1 •, 1 




ol^OJ — 


Lo/ 2 ;/ [7 1 c; _ lOiA-r- -i — ^Dt 7/1 ■> 1 


C A97r"l 


— 


2 

o/ 2 ;/ [Y 1 ^ ^D;/W •, -I- 'fdr 7/1 •, 1 

^ [I 10 ovvjx^k -+- oux<j«; : ,fc>J , 


I _r\Z ( (J. J 


c !TT('7'l - 
oii{t ) — 


L^/ 2 n _ O7/V — 1 "l-r- • ■> -1- ^Dt -i/i-,^ 

10 /V lox ijk + OUX <t Wjk) , 




V1(A) = 


1 T 2 (45 140z/ + 20z/ 2 )f 15x,-,- t + 30iwWiJ 


(A27f) 




— t 2 (2 — 9z/ 4- 2z/ 2 V — 15xvoi. 4- 30x^»-'Ltf*i.O 




TI(2) = 


^ 7 2 (8 - 21i/ + 2z/ 2 )(-15% fe + 30x <iWjk> ) . 


(A27h) 


Non-compact term: 








SI(5NC) = 2 7 2 (l + 2z/)% fc . 


(A28) 
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